{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# _*Pricing Short Butterfly Options*_ \n",
    "\n",
    "The latest version of this notebook is available on https://github.com/Qiskit/qiskit-iqx-tutorials.\n",
    "\n",
    "***\n",
    "### Contributors\n",
    "Stefan Woerner<sup>[1]</sup>, Daniel Egger<sup>[1]</sup>\n",
    "### Affliation\n",
    "- <sup>[1]</sup>IBMQ"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Introduction\n",
    "<br>\n",
    "Suppose a <a href=\"http://www.theoptionsguide.com/short-butterfly.aspx\">short butterfly option</a> with strike prices $K_1 < K_2 < K_3$, with $K_2 - K_1 = K_3 - K_2$ and an underlying asset whose spot price at maturity $S_T$ follows a given random distribution.\n",
    "The corresponding payoff function is defined as:\n",
    "<br>\n",
    "<br>\n",
    "$$ F(S_T) = \n",
    "\\begin{cases}\n",
    "0         ,& S_T < K_1 \\\\\n",
    "K_1 - S_T ,& K_1 \\leq S_T < K_2 \\\\\n",
    "S_T - 2K_2 + K_1 ,& K_2 \\leq S_T < K_3 \\\\\n",
    "0         ,& S_T \\geq K_3. \n",
    "\\end{cases}$$\n",
    "<br>\n",
    "In the following, a quantum algorithm based on amplitude estimation is used to estimate the expected payoff, i.e., the fair price before discounting, for the option:\n",
    "<br>\n",
    "<br>\n",
    "$$\\mathbb{E}\\left[ F(S_T) \\right].$$\n",
    "<br>\n",
    "The approximation of the objective function is explained in detail in the following paper:<br>\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "\n",
    "from qiskit import BasicAer\n",
    "from qiskit.aqua.algorithms import AmplitudeEstimation\n",
    "from qiskit.aqua.components.uncertainty_models import LogNormalDistribution\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariateProblem\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariatePiecewiseLinearObjective as PwlObjective"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Uncertainty Model\n",
    "\n",
    "We construct a circuit factory to load a log-normal random distribution into a quantum state.\n",
    "The distribution is truncated to a given interval $[low, high]$ and discretized using $2^n$ grid points, where $n$ denotes the number of qubits used.\n",
    "The unitary operator corresponding to the circuit factory implements the following: \n",
    "$$\\big|0\\rangle_{n} \\mapsto \\big|\\psi\\rangle_{n} = \\sum_{i=0}^{2^n-1} \\sqrt{p_i}\\big|i\\rangle_{n},$$\n",
    "where $p_i$ denote the probabilities corresponding to the truncated and discretized distribution and where $i$ is mapped to the right interval using the affine map:\n",
    "$$ \\{0, \\ldots, 2^n-1\\} \\ni i \\mapsto \\frac{high - low}{2^n - 1} * i + low \\in [low, high].$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# number of qubits to represent the uncertainty\n",
    "num_uncertainty_qubits = 3\n",
    "\n",
    "# parameters for considered random distribution\n",
    "S = 2.0 # initial spot price\n",
    "vol = 0.4 # volatility of 40%\n",
    "r = 0.05 # annual interest rate of 4%\n",
    "T = 40 / 365 # 40 days to maturity\n",
    "\n",
    "# resulting parameters for log-normal distribution\n",
    "mu = ((r - 0.5 * vol**2) * T + np.log(S))\n",
    "sigma = vol * np.sqrt(T)\n",
    "mean = np.exp(mu + sigma**2/2)\n",
    "variance = (np.exp(sigma**2) - 1) * np.exp(2*mu + sigma**2)\n",
    "stddev = np.sqrt(variance)\n",
    "\n",
    "# lowest and highest value considered for the spot price; in between, an equidistant discretization is considered.\n",
    "low  = np.maximum(0, mean - 3*stddev)\n",
    "high = mean + 3*stddev\n",
    "\n",
    "# construct circuit factory for uncertainty model\n",
    "uncertainty_model = LogNormalDistribution(num_uncertainty_qubits, mu=mu, sigma=sigma, low=low, high=high)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAEyCAYAAADOV2anAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3debgcVZnH8e9PEAgEMLKERSSACiJRMEGMoCSArDpBEKKiYxCNjAKOIgwiQgAdAWVxYBQjakRGoyKiIIhhSRDZ40JYomwhBhQEAjGEQELe+eOcC51K973dt/tWNcnv8zz93FtVp6re7tu33646myICMzOzgfaKqgMwM7OVgxOOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpKk84kraVdI2khZIekXSKpFX62OdNkn6Tyz8naY6kCyRtXCg3WVLUeWwzsM/KzMyKVq3y5JKGAFcDdwNjga2AM0mJ8IRedl0XeBC4EHgE2AI4CRghaceIWFJTdhZwaGH/2c3Et/7668ewYcOaKTognnnmGdZaa63Kzt9It8YFjq0/ujUucGz9UXVcM2bMeDwiNqi7MSIqewBfAOYB69SsOxZYWLuuyWO9GwjgrTXrJgO39ze+ESNGRJWuu+66Ss/fSLfGFeHY+qNb44pwbP1RdVy9feZWfUttH+CqiJhfs24KMAjYtcVjPZF/rtaJwMzMrLOqTjjbkG55vSgi5pCucPqsZ5H0CkmrSdoaOA24Dbi1UGxbSfNzXc8NklpNZGZm1gGKCsdSk7QYOCYizimsnwtcGBHH97H/b4C98uIMYN+IeKxm+2eA50l1RBsARwMjgF0iopiYevaZAEwAGDp06IgpU6b056l1xIIFCxg8eHBl52+kW+MCx9Yf3RoXOLb+qDquMWPGzIiIkXU3NrrXVsYDWAz8Z531c4H/bmL/1wM7AR8mXSnNANbopfyapMYGlzYTn+tw6uvWuCIcW390a1wRjq0/qo6LLq7DmUdqcVY0JG/rVUTcGxG3RMRFpCudHYAP9VJ+IXAF8Nb+hWtmZv1VdcKZRaGuRtJmpCuRWXX3aCAiHgKeBLbsq2h+mJlZiapOOFcCe0lau2bdOOBZYHorB8oNB9Yj3TJrVGYQsB/p1puZmZWo0o6fwPnAUcAlkk4nXZ1MBM6KmqbSku4DpkfEYXn568AS4BbgKeCNpP4795OaVSNpXeBy4CLgPmB94LPAJsBBJTw3MzOrUWnCiYh5knYHzgMuIyWPs0lJp9aqQO1wN7cDR5Jak60BzAF+Dnw1Ip7JZZ4D/kkasWBDYBFwE7BrRNw+EM/HzMwaq/oKh4i4G9itjzLDCstTyFcyveyzCDig3fjMAIYd9+u2j3H08CWMb/M4s0/br+04zKpSdR2OmZmtJJxwzMysFE44ZmZWCiccMzMrhROOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVJxxJ20q6RtJCSY9IOkXSKn3s8yZJv8nln5M0R9IFkjauU3aspJmSFkm6W9K4gXs2ZmbWyKpVnlzSEOBq4G5gLLAVcCYpEZ7Qy67rAg8CFwKPAFsAJwEjJO0YEUvy8XcBfg58EzgK2Bf4saR5EfHbAXlSZmZWV6UJBzgcGAQcEBHzgamS1gEmSjojr1tORNwI3FizapqkucBvgTcDf8jrvwRcHxFH5eXrJL0JODGXNTOzklR9S20f4KpCYplCSkK7tnisJ/LP1QAkrQ6MAX5aKDcFGCVp3dbDNTOz/qo64WwDzKpdERFzgIV5W68kvULSapK2Bk4DbgNuzZu3Al5ZPD5wD+l5v6G90M3MrBWKiOpOLi0GjomIcwrr5wIXRsTxfez/G2CvvDgD2DciHsvbdgZuAHaIiD/V7PM64F5gr3r1OJImABMAhg4dOmLKlCn9fXptW7BgAYMHD67s/I10a1wwcLHNfPjpto8xdBA8+mx7xxi+aecvzFfGv2cndGtsVcc1ZsyYGRExst62qutw2nUk8Grg9aRGBldK2jkiFvX3gBExCZgEMHLkyBg9enQn4uyXadOmUeX5G+nWuGDgYht/3K/bPsbRw5dw5sz2/uVmHzK67TiKVsa/Zyd0a2zdGhdUn3DmkVqcFQ3J23oVEffmX2+R9DtSy7UPAd+r2b94/CE15zYzs5JUXYczi0JdjaTNgDVZvu6lVxHxEPAksGVedT+wuHj8vLwU+Gs/4jUzs36qOuFcCewlae2adeOAZ4HprRwoNxxYj3SVQ0Q8B1wHHFQoOg64KSLavylvZmZNq/qW2vmkDpmXSDqddHUyETirtqm0pPuA6RFxWF7+OrAEuAV4CngjcCzpqqa2lv9UUh+dc4BLSR0/9wX2HtinZWZmRZUmnIiYJ2l34DzgMlLyOJuUdGqtCtQOd3M7qcHABGANYA5pRIGvRsQzNce/QdL7gS8D/0Gu4/EoA7aiGNahxgztNoqYfdp+bcdhK76qr3CIiLuB3fooM6ywPIVlr2R62/dS0tWNmZlVqOo6HDMzW0k44ZiZWSmccMzMrBROOGZmVgonHDMzK4UTjpmZlcIJx8zMSuGEY2ZmpXDCMTOzUjjhmJlZKZxwzMysFE44ZmZWCiccMzMrRcujRUsaDrwN2Ig0NcCTpNkzb4wIT9tsZmZ1NZVwJG1Jmk/mEGAoaYrmp4DngFeRpoReKmk6cAHwk4hYOiARm5nZy1Kft9QkXQDcBWwPnALsAKwRERtExGsiYjCwIfBeYCZwBnCPpF0GLmwzM3u5aeYK51lgm4h4qFGBiHgcuBK4UtLngIOATTsTopmZrQj6vMKJiCN7SzZ1yi+NiJ9ExE+aKS9pW0nXSFoo6RFJp0hapY99dpT0fUn35f3+IukkSWsUyk2UFHUeezf7fMzMrDPammJa0nbAroCA6RExs8X9hwBXA3cDY4GtgDNJifCEXnYdl8ueDtwLvBk4Nf88sFD2aaCYYO5pJU4zM2tfvxOOpP8AvgJcA6wFfE3S0RHxzRYOczgwCDggIuYDUyWtA0yUdEZeV89p+TZej2mSFgHflrR54YpsSUTc3EJMZmY2AJppNLBmg03/BYyKiIMiYl/gCOCLLZ5/H+CqQmKZQkpCuzbaqZBsevwx/9ykxRjMzKwEzXT8/KukQ+qsF6l5dI/+NIPeBphVuyIi5gAL87ZWjMox3F9Y/ypJj0taLOmPkg7oR5xmZtYmRUTvBaR3AecAzwNHRcStef2nSc2kryH1w9kdODYizm365NJi4JiIOKewfi5wYUQc3+RxNgLuAK6IiPE16z9MarL9R2Bt4JPAvsCBEXFJg2NNACYADB06dMSUKVOafTodt2DBAgYPHlzZ+Rvp1rhg4GKb+fDTbR9j6CB49Nn2jjF803WXWe7WuDplZXyvtavquMaMGTMjIkbW29ZnwgGQJOAwUsX8VOC/IuLvkt7CS7e+ro+IP7USWCcSjqTVSA0PXgOM6G20g/w8bgQGRcT2fR175MiRcfvtt/dVbMBMmzaN0aNHV3b+Rro1Lhi42IYd9+u2j3H08CWcObOtdjrMPm2/ZZa7Na5OWRnfa+2qOi5JDRNOU2OpRXIBsDXwKHCnpC8CsyLif/KjpWSTzQPqfTUakrf1KieQC4E3Afv2NbROpOx6CfDmvppem5lZZ7U0eGdEzI+IY4CdSOOpzZL0/jbOP4tCXY2kzUi36GbV3WNZ55CaU4+NiGbKA0R+mJlZiZpqpSbpy5JuyZXuk4BFETGWVNdxkqTp+fZaq64E9pK0ds26caTRDab3EdcXSC3jPhwRNzRzsnxFdCDw54h4oR/xmplZPzVz4/a7wLakPjcLSUlmqqRtI2JqTjSfyusujYgJLZz/fOAo4BJJpwNbAhOBs2qbSku6j9Sx9LC8/CHgv4HJwMOS3l5zzPsj4p+53HTg56SrpbWAT5CuzvZvIUYzM+uAZhLOPsBBETEVQNLvgSdIPf3vy6NCnyfpR8BJrZw8IuZJ2h04D7iMNAL12aSkU4yzts5lz/xzfH7UOpSUiADuA/4T2JjUZPoPwH4RcWUrcZqZWfuaSTizgI9ImgEsIjUtfgaYW1soIp4EPtNqABFxN7BbH2WGFZbHs3yiqbffYa3GY2ZmA6OZhPNR0hXD46TK9tmkK55FAxeWmZmtaPpMOBHxF2CUpLWA1Tyrp5mZ9UfTvb0i4hnSrTQzM7OWNdMs+iOtdpKU9DpJ7+x/WGZmtqJppuPn54D7JZ3aW18bSetJOkTSZcCfSC3DzMzMgObqcHaQNA44EviipAWkCcweB54DXgVsAbyWNBzNRcDhEfHwgEVtZmYvO03V4eTpon8iaStgD+CtwEakzpSPAtcDvwemRcTiAYrVzMxexloaIjYi7mf5+WbMzMz61NLgnWZmZv3lhGNmZqVwwjEzs1I44ZiZWSlaSjiS3ivJScrMzFrWavK4FJgr6XRJbxyIgMzMbMXUasLZCvgOcDBwp6SbJH1C0jqdD83MzFYkLSWciJgdESdFxBbAu0kTnJ0N/F3SDyWNGYggzczs5a/f9TERcW1EfAR4AzADOAS4WtIDkj4rqaVOpWZmtmLrd8KRtKukycBfgO2A/yVN/XwxcDJwYScCNDOzFUOrrdQ2l3SipPuBa4HNgAnAxhFxZERcExHHkmYJHdvkMbeVdI2khZIekXRKX9MhSNpR0vcl3Zf3+4ukkyStUafszpJukbRI0oOSjmrlOZuZWWe0etvrAeAR0pTT34uIBxuUuwu4ta+DSRoCXA3cTUpQWwFnkhLhCb3sOi6XPR24F3gzcGr+eWDN8V8HXAVcDnwBeBtwlqSFEXFBX/GZmVnntJpw3gNcFRFLeysUEX8FmmlAcDgwCDggIuYDU3OLt4mSzsjr6jktIh6vWZ4maRHwbUmbR8RDef0xpAT54YhYAlwr6bXASZK+GxHRRIxmZtYBrdbh7EialmA5kjaWdGKLx9uHlMBqE8sUUhLatdFOhWTT44/55yaF41+Sk03t8V9DqncyM7OStJpwTiJ9WNezSd7eim2AWbUrImIOsDBva8UoYCl5+gRJa5HqmGYVyt1Tc24zMyuJWrmrJGkpsFNE3FZn21jguxGxfgvHWwwcExHnFNbPBS6MiOObPM5GwB3AFRExPq/bFJgLvC8iLq0puyqwGPhkREyqc6wJpIYQDB06dMSUKVOafTodt2DBAgYPHlzZ+Rvp1rhg4GKb+fDTbR9j6CB49Nn2jjF803WXWe7WuDplZXyvtavquMaMGTMjIkbW29ZnHY6kj5JanQEE8C1JxbqVNYDhwG/bCbQ/JK0G/BRYAHy23ePlJDQJYOTIkTF69Oh2D9lv06ZNo8rzN9KtccHAxTb+uF+3fYyjhy/hzJntdU+bfcjoZZa7Na5OWRnfa+3q1riguUYDC4En8u8CngaeLJR5HrgS+GaL558H1PtqNCRv65Ukkfr7vAnYOSJq93kq/ywef0jNuc3MrCR9JpyI+BnwMwBJ3wdO6aU5dKtmUahLkbQZsCbL173Ucw6pOfW7I6JYF/SMpL8Vj1+z3MzxzcysQ1odS+3QDiYbSFdFe0lau2bdOOBZYHpvO0r6AnAEqcnzDb0c/32FjqTjgL8Bd/Y7ajMza1nVc9ucDzwHXCJpj1xhPxE4q7apdB5R4Ls1yx8C/pt0O+1hSW+veWxQc/yvkVrV/VDSGEnHAp8kXaW5D46ZWYmaaTRwKzA+Iu6WdBup4UBDEfG2Zk8eEfMk7Q6cB1xGqnc5m5R0inHWXqXsmX+Oz49ah5JGQiAi7pO0N3AW6WrnH8DRHmXAzKx8zTQauIt0i6vn945eGUTE3cBufZQZVlgez/KJptG+N5CGtDEzswo102jg0Jrfxw9oNGZmtsKqug7HzMxWEs3U4fRZb1OrlTocMzNbeTRbh+MWXWZm1pZm6nDGlxCHmZmt4FyHY2Zmpai0H46Zma08Ku+HY2ZmKwf3wzEzs1K0PAlGnn9mPKn3/sbA34FbgB9ExPMdjc7MzFYYLTUakPRG4F7gf4HtgBfyz/8F7pO0bccjNDOzFUKrVziTSBOwvTMi5vSslPRa4HLS6M/v6lx4Zma2omg14YwEPlibbAAiYo6kk4AfdSwyW+kM69B0ye1Ouzz7tP3ajsPMltdqP5zZwBoNtq0BzGmwzczMVnKtJpzjgC9L2ql2paS3A6cC/9WpwMzMbMXSn8E71wFulPQY8BiwYX48ARwPXDoAcZqZ2ctcfwbvvGuAYjEzsxVY5YN35qbU5wKjSFNMXwCcHBEv9LLPasBXgLeTGjKsERGqU24y8NE6h3hjRMxqP3ozM2tWyx0/O0nSEOBq4G5gLLAVcCapbumEXnZdE/g4cCtwI71PUT0LOLSwbnb/IjYzs/6qNOEAhwODgAMiYj4wVdI6wERJZ+R1y4mIpyS9OiJC0hH0nnCeiYibOx+6mZm1ouXpCSSNk3S1pDmSHis+WjzcPsBVhcQyhZSEdu1tx4jwIKJmZi8jrQ5t8yHgB8B9wGuAX5FGGHgFMB84r8Xzb0O65fWi3Kl0Yd7WCdtKmi/pOUk3SOo1kZmZ2cBQKxcKkv4IXAycBiwGRkbEHyStDUwFLo6Ir7dwvMXAMRFxTmH9XODCiDi+iWMcAZzboNHAZ4DnSXVEGwBHAyOAXSLi1gbHmwBMABg6dOiIKVOmNPt0Om7BggUMHjy4svM3MlBxzXz46baPMXQQPPps3+V6M3zTdZdb162xdWtcndKt/wPQvbFVHdeYMWNmRMTIettarcN5PfD7iHhB0gukPjlExL8knQ6cDTSdcAZaRHyjdlnSFaRm3ccD+zfYZxJpzDhGjhwZo0ePHuAoG5s2bRpVnr+RgYqr3SFpIA1tc+bM9qomZx8yerl13Rpbt8bVKd36PwDdG1u3xgWt1+HMB1bPvz8MvLFmm4D1WjzePKDeV6MheVtHRcRC4ArgrZ0+tpmZ9a7VrzW3AW8GriLV35woaQnpttWJQKutwWZRqKuRtBmp2fNA9ZMJPGupmVnpWk04XwU2z7+fmH//FulK6Tbgky0e70rgGElrR8S/8rpxpCmtp7d4rD5JGgTsB8zo9LHNzKx3LSWc3J/l5vz7U8BYSasDqzfqM9OH84GjgEtyHdCWwETgrNrjSboPmB4Rh9Ws2wdYC9g+L78/b7otIh6StC6pBd1FpFZ16wOfBTYBDupHrGZm1oaOTTEtqeUppiNinqTdSc2pLyMNbXM2KekU41ylsO5bvHS1BfCz/PNQYDLwHPBP0ogFGwKLgJuAXSPi9lbiNDOz9rWUcPIU078hXSXMII0WvR3w78CXJO0dEXe3csxcvreRAoiIYc2sK2xfBBzQSixmZjZwPMW0mZmVotVm0SOBE+tNMQ2cBOzYqcDMzGzF4immzcysFK3eUjsOOFPSgxFxS8/KmimmP9/J4Mzs5WtYh0ZBaHc0hdmn7dd2HNYZnmLazMxK4SmmzcysFJVPMW1mZiuHfg0RK2kTYBTwatKttJsj4pFOBmZmZiuWVjt+rgKcC3yCZXv+vyBpEnBkRCztYHxmZraCaLVZ9MnAx0iNA4aRpoIelpc/xvJD0piZmQGt31L7d+CEwqyec4CvSQrSQJwndio4MzNbcbR6hbMhcEeDbXfk7WZmZstpNeH8FfhAg20fAP7SXjhmZraiavWW2peBKXmwzouBR0lXNQcBY2icjMzMbCXX6gRsP5X0FKnxwDeAVwKLSVMV7B0RUzsfopmZrQiaTjiSXkmadO3OiBgl6RWkWTQfd1NoMzPrSyt1OC8A1wLbAETE0oh4zMnGzMya0XTCyYnlXmCjgQvHzMxWVK22UvsicKKk4Z0KQNK2kq6RtFDSI5JOySMa9LbPapK+Jul3kp7NfYAalR0raaakRZLuljSuU7GbmVnzWm2ldgKwHvAnSQ+TWqkt82EfEW9r9mCShgBXA3cDY4GtgDNJifCEXnZdE/g4cCtwI7Bbg+PvAvwc+CapU+q+wI8lzYuI3zYbp5mZta/VhHMXcGcHz384aXicAyJiPjBV0jrAREln5HXLiYinJL06IkLSETRIOMCXgOsj4qi8fJ2kN5FGQ3DCMTMrUavNosd3+Pz7AFcVEssU4HRgV+CyXmJpeBsNQNLqpL5BRxU2TQG+L2ndiHi6X1GbmVnLmqrDkTRI0oGSjpb0IUlDO3T+bYBZtSsiYg6wMG9rx1akfkKzCuvvIT3vN7R5fDMza4H6uFBA0pakepZhNavnAwe3Ww8iaTFwTEScU1g/F7gwIo5v4hhHAOdGhArrdwZuAHaIiD/VrH8dqbXdXvXilzQBmAAwdOjQEVOmTGn9iXXIggULGDx4cGXnb2Sg4pr5cPsXnEMHwaPPtneM4Zuuu9y6bo2tW+OC7o6tE1a2/89mjRkzZkZEjKy3rZlbamcAS4F3kkYU2IJUCf/t/PsKJSImAZMARo4cGaNHj64slmnTplHl+RsZqLjGH/frto9x9PAlnDmzX/MKvmj2IaOXW9etsXVrXNDdsXXCyvb/2QnN3FIbRZqS4PcRsSgi7gE+CbxW0sZtnn8eUO/rx5C8rd1jU+f4QwrbzcysBM0knI2BBwrr7gdE+51AZ1Goq5G0GanZc7HupVX3k8Z5K9YFbUO6Yvtrm8c3M7MWNNvxs/eKnv67EthL0to168YBzwLT2zlwRDwHXEcaybrWOOAmt1AzMytXszdHr5K0pM76a4rrI6KVSdjOJzVbvkTS6cCWpGmqz6ptKi3pPmB6RBxWs24fYC1g+7z8/rzptoh4KP9+KjBN0jnApaSOn/sCe7cQo5mZdUAzCefkgTp5RMyTtDtwHqnPzVPA2aSkU2tVoDjczbeAzWuWf5Z/HgpMzse/ISeiLwP/ATwIfMijDJiZla/PhBMRA5Zw8vHvpvFIAT1lhjWzrsG+l5KubszMrEKtDt5pZmbWL044ZmZWCiccMzMrhROOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVJxxJ20q6RtJCSY9IOkVScTrpevutK+n7kuZJelrS/0lar1BmsqSo89hm4J6RmZnV0+cU0wNJ0hDgauBuYCywFXAmKRGe0MfuPwXeAHwcWAqcTppK+p2FcrOAQwvrZrcTt5mZta7ShAMcDgwCDoiI+cBUSesAEyWdkdctR9IoYE9g14i4Pq97GLhF0h4RcXVN8Wci4uaBfRpmZtaXqm+p7QNcVUgsU0hJaNc+9nu0J9kARMStwIN5m5mZdZmqE842pFteL4qIOcDCvK3p/bJ76uy3raT5kp6TdIOk3hKZmZkNEEVEdSeXFgPHRMQ5hfVzgQsj4vgG+00l3Srbv7D+ImDLiHhHXv4M8DypjmgD4GhgBLBLviKqd+wJwASAoUOHjpgyZUobz7A9CxYsYPDgwZWdv5GBimvmw0+3fYyhg+DRZ9s7xvBN111uXbfG1q1xQXfH1gkr2/9ns8aMGTMjIkbW21Z1Hc6Aiohv1C5LugK4Czge2L/BPpOASQAjR46M0aNHD3CUjU2bNo0qz9/IQMU1/rhft32Mo4cv4cyZ7b2tZx8yerl13Rpbt8YF3R1bJ6xs/5+dUPUttXlAva8fQ/K2ju4XEQuBK4C3thCjmZl1QNUJZxaFOhdJmwFrUr+OpuF+WaO6nVqRH2ZmVqKqE86VwF6S1q5ZNw54Fpjex34bSdqlZ4WkkcCWeVtdkgYB+wEz2gnazMxaV3XCOR94DrhE0h65wn4icFZtU2lJ90n6bs9yRNwE/Ba4UNIBkvYH/g+4oacPTh6J4HeSPilpd0njgOuATYD/LusJmplZUmmjgYiYJ2l34DzgMuAp4GxS0qm1KlAc7mZcLvs9UuK8HDiqZvtzwD9JIxZsCCwCbiJ1Fr29o0/EzMz6VHkrtYi4G9itjzLD6qx7ijRkTXHYmp7ti4ADOhCima1ghnWoBV07LfFmn7Zf2zG83FR9S83MzFYSTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVjzRg5eqGHtawcvayNlvZ+QrHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZWi8o6fkrYFzgVGAU8BFwAnR8QLfey3LnAOsD8pcV4OHBURTxTKjQW+DLweeCAf+yedfh5mZu1a0TtmV3qFI2kIcDUQwFjgFOBo4OQmdv8pMBr4ODAe2BG4tHD8XYCfA9cB+wC/Bn4sac+OPAEzM2ta1Vc4hwODgAMiYj4wVdI6wERJZ+R1y5E0CtgT2DUirs/rHgZukbRHRFydi34JuD4ijsrL10l6E3Ai8NuBe1pmZlZUdR3OPsBVhcQyhZSEdu1jv0d7kg1ARNwKPJi3IWl1YAzpSqjWFGBUviVnZmYlqTrhbAPMql0REXOAhXlb0/tl99TstxXwyjrl7iE97zf0I14zM+snRUR1J5cWA8dExDmF9XOBCyPi+Ab7TQWeiYj9C+svAraMiHdI2hm4AdghIv5UU+Z1wL3AXhGx3G01SROACXlxa+Av/X6C7VsfeLzC8zfSrXGBY+uPbo0LHFt/VB3X5hGxQb0NVdfhdJ2ImARMqjoOAEm3R8TIquMo6ta4wLH1R7fGBY6tP7o1Lqj+lto8oF5dypC8rZ39en4Wyw0pbDczsxJUnXBmUairkbQZsCb162ga7pfV1u3cDyyuU24bYCnw137Ea2Zm/VR1wrkS2EvS2jXrxgHPAtP72G+j3M8GAEkjgS3zNiLiOVL/m4MK+44DboqIp9sPf8B1xa29Oro1LnBs/dGtcYFj649ujavyRgNDgLuBO4HTSQnjLOCciDihptx9wPSIOKxm3VWk0QM+T7piOR14LCLeWVNmF2AacB6pU+i+ufze9RoMmJnZwKn0Cici5gG7A6sAl5FGGDgbOKlQdNVcptY40lXQ94ALgRnA+wrHvwF4P7AHcBXwb8CHnGzMzMpX6RWOmZmtPKquwzEzs5WEE1I1XdsAABkPSURBVI6ZmZXCCcfMzErhkQa6gCSRGjzsB7wReDWp5d0/gJuByRFRSb+h3C9qX0DAzyLiCUmvIbX22wqYDUyKiJklxvRfwBVlnrNZkgYBq0bEv2rWbQAcAWxL+rv+Cfjmy6RpfiXy/8R7gbeSpi+5nfQ374pK5zyq/ePAbrlxUlUx7AasBvw6Ip7J77VPk1r8PkD633ykivjqcaOBiuU3yBXACFKCeR7YlPRPdiXpjbM1cGpEnFpybG8DpgJrAUuAJ4G9crwvAHcB2wEbAXtExO9Kimsp6fWZBfwImBIR95dx7r5IugK4NyI+k5dHkf6OS0ktKUX6Wz9P+rC6q6S4dgAGRcSNNev2Br7AS4nwz8DE2jIlxXYjcFhE3JOXh5CmDxkBLMjFBpO+fO1Vm8wHOK5P9bJ5EPA14BuksRmJiG+WERe8OCbkNcBmedWDpClbpgKvInV835rUp3FERMwtK7ZeRYQfFT6AH5PesMNr1m0C/Ab4eV7elfSP97GSY5tK6jz7KtLI2+cBc4FfAq/MZVYnfaBeV2JcS4GvkmZ5fY6U/G4DPgtsWvHf83FgbM3yzaQPhrVr1q1LatJ/VYlx3Qx8sWb5Y/l1vAb4InBC/lsvqY2/xL/n22qWv0v6crN3zbq9ScNRnV1yXC/kn/UetdteKPk1+ynpC8LrSHdEfpg/R27sea+RBvH8M/DtMmPrNe6qA1jZH6RptQ+ss35YfkNvnJePB/5ccmxPAPvULG+Y/7n2LJTbD3i8xLhe/IDK/2wT8gfnkvyYlte9uoK/50LgXTXLzxdfr5rX7JkS45pfGwdwH3BunXLnV/A+KyacfwL/Wafc54GHSozrF8DfgUPJd4Nqtr0qx/2usuIpnP8R4OCa5c1zPAcUyh0K/LWKGOs93Gigeq8gJZaiF0i3X3oGH72FaubwiTq/F+/DVnZfNiKejIhJEbE78BrSFOWrkT44/y6p/UniW3MnaeK/Ho+SkmLReqTkVJalheXNgYvrlLuYdCumSq8i1dkUzSDdvi1FRLwP+ChwDHBbnvLkxc1lxdHAENIt+B4P558PFco9QPq/6ApOONWbCnxZ0pY9K/I97P8hvaF6GgsMBsquZL4d+LykwZJeQbrKehj4D0mr5FhXBT5F+qCtVET8IyK+ERHvALYgjVixSclhnAYcJ+lj+bX5CvA1Se+WtJqk1XPdyVdZfjbagfQ74JCa5buAekPY78hLH15lOlDSp3K9yTyg3nwq65Ou1EoTaVSSN5NmCv61pCm50UzVHiN9aejxAvBt0hecWhsCpdR5NaXqS6yV/UH69nEXaWTr+0hjyz1LutVWezvrDOAnJcc2kvTPvzjH9ATwFlISfIA0HNGDpHqUMSXGtcwtmG57AB8nfTA+Ddyaf3+BdLvvhfz4BbBmiTENz3H8EHgbaSr2x0gJ8d2kCufTgEXUuZ1Vwt+z+PhenXLfBn5X4d91I9IwWguAM/PfsapbapfWe43qlDsXmFrVa1Z8uJVaF8hXCweTPszXICWeH0XEk5UGBuRvc+8hNaH/eUT8XdJGwLGkWy8PARdExB9KjOkk4DvRRc09iyStRxrv722kDyqRkvc9wOURMaOCmLYHvgXsRLolpLyp5/d5wCkR8Y2yY2uGpE8A90fEtRXHMYo05uPWwH5Rcqu+HMNQ0heWB/so9zlSndw15UTWOyccs5WMpDeSkk4xEd4YEYurjM1WbE44XUTSm0gTxNXOSjorSuqr0SpJr4iIYmV0ZSStQeqMuhS4r+oPz1yHsyU1HXkjYk6VMb3c5A6gRIUfVLkzryJiYc267ckdn6u4Wn25cqOBLpArmB8C7gB+RppAaVL+/Q5JsyUdWlFsB0i6VNIVkt6b142TNBtYLOmhfKujzJg+LOljNcurSjqN1HfjDlIDhiclHVdmXDXxjJD0K1Jl7T3A74GbgAclPSzpFElrVhFbN5K0Z2ESRiTtL+kPpPrD5yXdLmm/kuNaV9IvSHVf8yV9R9Iqkn4A/IH0/3mrpBskrV9mbM2SdKCkeq1gK+GEUzFJR5IqQy8HRpNalbwyPzYkdfq8HDhf0qdLju1gUjPZ9Un/+D/JyeWHpH4vR5E6mp0vaa8SQzue1OG0x+k5lq8C7yK9ZmcCJ0k6vsS4kLQn6TXZhHSf/1RSr/kXgImkCQYPBG7MrRHLjO09kq6RdI+kX0p6V50yO1XwAXUlaUinnhjeB1xCasBwXH4sBn6ZX9+ynAq8E/gcqaPsO0gtC3cjdUQdSqrf3CKXtT74llrFJD0AnB8RZ/RR7ljg8IjYsrdynSTpNmBGRByelw8hTXh3XkQcXVPu+8BmEbFHSXEtJLXgm56XHwO+UqzslvR54MiI2LzOYQYqthnAnRHx0cL6I0l9hLYk9RO6Ebg5InobPqWTcb2bNHrFzcAfgVHA9sA5wOd7bllJ2olUl1Oc8HAgY1sKvD0ibs3LfwAejoj3FspdAawVEbuWFNds0vvqO3l5B1JfoEMj4gc15T4BHB8RW5QRVz7n95osujkwusy/Z298hVO9jUlNZ/tyKyV2esu2ZtnOgZeTrryKnSkvIY3HVZanSVddPdYlDeFR9GfSVWKZtgUuqrP+IuC1wNYRsYj0Qf++OuUGyknAhRGxc0QcEREjgE8AnwQuyfVf3WI70lV/0STSYJ5l2YCX+sFBHjONNE5Zrfuo329oIH2U1JR9eB+P0r5sNcMJp3p/Bj6RO1bWlStOP0GqnyhTsOzU3j0DKT5VKLeA1Du8LL8idUhdLS9fDXywTrkPkvo4lekxUvP2oreQXs+ezrsP8dIoEmXYjkIijIjvkW4/vh24VlK9ERHKUnur5Wnqd1Z8hnI/sx4kvT493klq/PGOQrmdgbIbg9wLXBsRO/b2IN2O7BqenqB6R5Nuddwt6RLSCMg9H+jrklqtvY/UQXTvkmN7iPSN/SqAiHgh90G4p1BuK9KYU2X5Aqnn/J2SLiB1QD1d0nakcdREGl5mB9IQ92WaBJwqaTApET5P6r3/RdIApz19h7ak3A+pRaRRv5cRETPykC1XkW7zTSwxplpXSVqSf1+X9LebXiizDeW+z84HviFpOCkJHkx6730p/33/TLri+izl1+HczPKJr57a/laVc8KpWET8PjexPJY09MhmhSJ/I1Wqfi3KH4L/EgrjMEXELXXKfQAobU6QiHhS0ttJH+Kf46XbZqPy43nSkEHvjIjbyoorx/aVXCdxHHBiz2rSqOD/WVN0MfDfJYZ2B2l0gV8VN0TEAznpXAFMLjGmHifXWfdYnXUHkka0LkVEnJfvPHyQ1DDg2Ig4X9Jc0tBTPePhnQ98vay4snNJLeX6Mp1lx/arlBsNdJncXLbn9tRTtW3/u5Wk15JiLXWcq5rzD2PZToz3d0EfnFeSrvzWAB6o6rWpieeTpNZ9OzQawULSWqQhd/aICN9u70W+zb1+RPyz6lheTpxwzMysFL6l1iWUpnLehPTt/PE629cH9o2IC0sPro58D/sPwCFl37ZSl0/jrC6clrvb6WUyXXK+sqmd+noGKd7Sv7lLGkm6zSjSNPSzJL2FdIuy5312XkRcVXZsjfgKp2KSVie1Hjogr1pKGpH2c7UflhX1j9i3l81rAT8h1VXcCRARV5QUV1dO45xj6cppuZulNM7aQRFxSonn7MrpktWlU1/nWPYiNZZ5ktR6bwNgLKne9R5SX6sRpAYrB0bEpWXF1qsyhqT2o9fhw08ktUr7BGk6gKNIc1rcC7y+ptxOlD+NbVdOsUuXTuOcz9uV03K3EP+BFbzPunK6ZLp06ut83t+ThtZZJS8fn+P4bqHcD0kdjCt/b0V4iunKH6Rm0EcU1m0EXE+aandUXldFwrmdl6bY3bzweHP+hzy4Z12JcXXlNM75nN06Lfdrm3wcXsH7rCunS66TcLpi6ut8zqdJV8g9y0NyvLsVyu1JatBTWmy9PVyHU73NKHTojIh/SNqd9O3k6jykTJn9D3rsSLryOp3U7+XzkeffkNTTafEfEVGc1nag9UzjfH1e7pZpnHt047Tcs5s8p5os10kvl+mSu2Lq6+xZlu1X1fP7oEK5NUl9sLqCE071HgFez0sfngBEatb7AUnfIF06l95YINJXpEmSfgp8GZgp6dz8e5VOA/5P0t9Ir0vPNM5PkG6j9XT8LHsaZ3hpWu4bSMmudlruayN1nq1iWu5/AdcCF/RRbhdSn7AydfN0yQfmynnooqmvSbfUTpR0bz7310mzBf+XpOsj4l/5S+GxpITYHaq+xFrZH6TBMKf1UeYLlFxP0iCON5N68j8CfIZqp9jtummcc1zdOi33VNJQKH2Vq6IOpyunS6aLp74m1XfNrnmv30+6JdrzvzCTlJznAduXGVtvD1/hVO+bwDhJr44GHfIi4qtK8+W8u9zQlovjDmC0pA8AZ1DhkBkRcUGeq6RnGucn6YJpnCPi9jwUSnFa7nfx0rTcV1LytNykK+gJTZT7J4Wr7RJ8knTrpy8PkpJTKaL5zq+3k1psliYi7stDOe1MapxyTUQ8K2k06cvY1qRb8j+Kklr1NcPNoq1f8m2htYAFEdE1EzyZWfdywjEzs1J4vKSXiTy97XerjqOebo2tW+OC7o6tW0m6WtI1VcdR1K1xQffF5jqcl48xdO8XhG6NrVvjgi6NTdLVpDsfu1cdSx2iC18zujcu6LLYfEvNzF6Uvw2/IiK6Zkh7W3F0Teaz3klaI08D0HW6NbZujQu6N7aI2L1bk42kV3bja9atcUH3xeaE8/KxH6lZaDfq1ti6NS7o0tiq+oCS9GlJ90t6VtKfJX2kTrG3UvJr1q1xdXtsjTjhmK0kuvUDKvfrOpc0COuXSJ0YJ0u6WNIaZcbycoir22PrjRsNVEzStU0WrTekxoDq1ti6NS7o3thqPqB+TBq6/h2kD6ixwIcjosrxtj4PfD0iXhxSJ48l+H/AdZLeExFPOK6XTWwNudFAxSQtAf5CGgepN5sCO0W58+F0ZWzdGhd0b2ySbicNbVPvA+pB4D2RJoqrYt6lfwHvjYhphfXDSKMyrEKaBmCDMmPr1ri6Pbbe+AqnencBsyJiXG+FJL2fkofPoHtj69a4oHtj25r0rfhFEXGNpLeTPqBukrR3ifHUepo0AOYyImK2pHcAvwZuAk51XC/q5tgach1O9W4G3t5EuaD8scu6NbZujQu6N7aGH1Ck22uPkz6gdiwxph4zgP3rbYiIecDupPHK/qfMoOjeuKC7Y2vICad6ZwBHNlHuCmCLAY6lqFtj69a4oHtj6+YPqIuALSXVm9OIiHgW+DfS1ApzHBfQ3bE15Docs5WApIOAz5LqauqOSi5pFeBbwLsjouxEbSsBJxwzMyuFb6mZmVkpnHDMzKwUTjhmZlYKJxwzMyuFE471StJ4STMk/UvSPEl/lHTWAJ3rYEnjmyg3UVLUPB6R9HNJWzV5nsm5533lmn3OuWzP8763wfZ78/aJAxVDi8dd5nXu9HkkvULSEfk9+ayk+ZLukvQ/kvrVx0nJnyR9tMH2ybk3f71t58mT6vXKCccakvQFUjv+q4ADgH8Hfklq3z8QDgbGN1n2aWBUfnwe2B64RtJaTex7agvnGWitPGeARcAWkkbWrpS0IzAsbx/oGJpVfJ07fZ6fAF8GLiG9Jz9K6t/0juh/89uDgVcDP+rHvl8HDpH0un6ee4XnoW2sN0cA346I42vWXSbp5KoCqrEkIm7Ov98saQ7wO2Bf4GfFwrmPySoR8XxE3F9inJ32DPAH4AOkjpo9PgBcC4yoIqgeZb3OkvYB3g/sGxFX1mz6RX+vbrKjgB9GxOKac61KSp4fATYBPijpfuDkiHhxeKI8rMwNwH8AR7cRwwrLVzjWm1cB/yiurP322HPbRNL+kmZJWiTpBknbFvfLt1RmSnpO0t8kfSX/MyNpMnAgsGvNrbKJLcQ6I/8cVieuu0jf/Heq3VaI7V2SrpO0QNLTkqZJ2qFm+zslTZe0UNITkr4jae3eApI0StKvJP1d0jP5Vs0hta9dP5/zFODgng/W/PPgvL5jMeTX4OLC8UbnMtvVvpZ9vc6NziNpX0lLJW1ROM8Wef3YBq/BrvnncqNz9/fqJl+ZvAO4uLDpM8CxpFEYrgA+BnwPWK/OYX5OusrxZ2sdvsKx3vwBODJfPVzey3DnmwNnkebleBY4GbhK0ut7hr2XtCfpFsiFwDHAm0nfGtcDDs+/v5aU5D6Vjzu3hViH5Z//KKw7Azglr687z4uk0cBU4DrSbZlngJ1JIzr/UdLOwNXApaRv1esBpwFD8nIjmwO/B84nfRDvDHxf0tKI+DH9f86XkEYE2IV0VfdO0qjAlwBfKymGWsPo+3VudJ6/A4+QXveJNeXHA4+RBqGs55n882uSzoyIh1qMuZ7d83H/XFi/K2mk7TPyF6nf5zHo6rkRGAoMr3Mciwg//Kj7ICWFB0gDTS4ljYR8CrBOTZnJefs7atZtDiwBDq9ZdzNwXeH4xwIvAK/JyxcD05qIayJpsMlV8+MNpGQxH9i4ENf2dfafDNxes3wT6faUGpzvd3Vi3y0ff7smX0vlWL9N+vDqWd/Uc6593vn3XwL/m3//JnBp/v1xYGInYgCmARcX1o2ufd4tvs6NzvNlUpJSTZyzSfO9NHotNgLuyOcO4E7geGBwG+/3ScBtddZ/G/hbPudkYFgvx1g1v/c/0d84VuSHL/usoYi4A3gjqUL2m6QPgi8Bt0saXFP0sYi4sWa/h0i3uN4GL97XfyvL1638hHRbd1Q/wlsPWJwffwG2BMZFxN9ryjwcEX/q7SC5kcFOwA8if2IUtq+Z4/uppFV7HsAN+dwN60wkDVFqMfVQTawTSAmyXVOA90tanXSVtdzttBJi6NHn69yH75G+pIzOy2Py8vcb7RAR/wB2APYiXe29CvgKcKOk1QAkfTTfQvxTvo07K/8+Q9Ir6xx2I1LCLvoK6crnQdL/wufzVW+9uJYAT+VjWYETjvUqIp6LiMsi4oiI2Bb4OPB64LCaYo/V2fUxYOP8+/rAK4FHC2V6luuOeNuHp0lD6Y8EXkP61nlloUzxfPUMISXSv/eyfRVSwl1c83iO9Jw26+XYk4FxpNtce+Z4vwd0YgrgXwGDSR+GawGXVRBDj2Ze54Yi4gHS1dShedWhwK0RcVcf+70QEb+NiE+Rbtd9n3Qra1Te/oOI2J70ZWcJsHNEbB8RI6KmUUCNNUh/1+J55uTjvo90xb8LcIMadw94js6+visM1+FYSyLiu5LOALapWb1hnaIbkm7BQfrWuLhOuaH5Z93Ri/uwJCL66kvTTOXxPNLtwo0bbH8qH2ciqcK46JF6OynNK/8e4NMRcX7N+o58yYuIZyRdThoB+mcR8UyxTAdiWASsVlg3pF44TR6vNxcA31Fqin8ALbbyioilkn5LSlbFD/vXA/Oi7ymXn6TBlUlOUL9Rmqp7Immqh7MlnZMTUq1X0b/39ArPVzjWkKTlEomkDYB1WfZb7YZKswz2lHkt6VvlrZC+iZJusR1UONzBpA/7m/Ly85T8zTB/UN8C/HtPq686228Gto6I2+s86iYcYHXS/9eL35hzq7ZiH6Z2nvO3SFc25zfY3m4Mc1n2iwWkq6T+6u25XpK3TyHFXPcWIYCkoQ02/RuwkPT3rPUWmqvA/wt15iiq974Abss/X10ouwGwJvDXJs630vEVjvVmpqRfAr8l3SLbnNTJciHwg5pyjwMXSTqBl1qpPUa6ndPjJFLLte+TPkyGk1oufScielpFzQLGStqf9GH3SC8f6J10HKkV2pWSJpHu148iVXhfTmrccI2kpaSK73+RbuHsB3wxIpb7cImIpyXdBpwoaT4psR5HuhW4Tk3Rfj/nSPPZT+tle7sx/AI4TNLZpNZiY4B2pqFu+FwjYpGk/wM+Dfw4Ip7q5Tg/lfQv4KekxgUbAocAY0mV9cV930JqYNCX35Neqw0i4p81638k6Y/A9aTblyNIV5YPA/cUjjGSdMV3I7a8qlst+NG9D9I//29Jt40Wkf65fwRsU1NmMqmF1wGkb3XPkf5xl2u9RapLmEn6JjuXVP+was329Ukfck+Sb2M1iGsiubVWL7FPpqaFVF/bSE1frycl06dIrd62r9m+E/AbUku4Z4C7SU3B1+0lhtcB1+Tyc0iJa5nYm33OLTzvZVqptRsD8AVSC61/kWaZ/DeWb6XW1Ovc13MF9sjr9+jjOX4s/y3m5vfSk6SEOLpB+cuADzTxfl8NeAL4SGH9+/L5/kFK2vNJiX6HOsf4BoUWjX689PAEbNaW3KFvu4gY2VdZs97kusGDgS0jYmkHjzsH2Csiilcj9cp+A3hdROzXYPtkUqKcXWfbKsBDwHERcVFbQa+gfEvNzColaWtgW9KQMCd3ONkMIXWKbbZO5WvAXyW9IercKu3DQaRbyg3rn1Z2bjRgZlX7NulW7RWk4WM6JiLmRcSgSA1Xmik/l3TLrlGrxUtJt1zrEXBYpL44VodvqZmZWSl8hWNmZqVwwjEzs1I44ZiZWSmccMzMrBROOGZmVgonHDMzK4UTjpmZleL/AW46+1xJ9j1lAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot probability distribution\n",
    "x = uncertainty_model.values\n",
    "y = uncertainty_model.probabilities\n",
    "plt.bar(x, y, width=0.2)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.grid()\n",
    "plt.xlabel('Spot Price at Maturity $S_T$ (\\$)', size=15)\n",
    "plt.ylabel('Probability ($\\%$)', size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Payoff Function\n",
    "\n",
    "The payoff function equals zero as long as the spot price at maturity $S_T$ is less than the strike price $K_1$, then increases linearly, and is bounded by $K_2$.\n",
    "The implementation uses two comparators, that flip an ancilla qubit each from $\\big|0\\rangle$ to $\\big|1\\rangle$ if $S_T \\geq K_1$ and $S_T \\leq K_2, and these ancillas are used to control the linear part of the payoff function.\n",
    "\n",
    "The linear part itself is then approximated as follows.\n",
    "We exploit the fact that $\\sin^2(y + \\pi/4) \\approx y + 1/2$ for small $|y|$.\n",
    "Thus, for a given approximation scaling factor $c_{approx} \\in [0, 1]$ and $x \\in [0, 1]$ we consider\n",
    "$$ \\sin^2( \\pi/2 * c_{approx} * ( x - 1/2 ) + \\pi/4) \\approx \\pi/2 * c_{approx} * ( x - 1/2 ) + 1/2 $$ for small $c_{approx}$.\n",
    "\n",
    "We can easily construct an operator that acts as \n",
    "$$\\big|x\\rangle \\big|0\\rangle \\mapsto \\big|x\\rangle \\left( \\cos(a*x+b) \\big|0\\rangle + \\sin(a*x+b) \\big|1\\rangle \\right),$$\n",
    "using controlled Y-rotations.\n",
    "\n",
    "Eventually, we are interested in the probability of measuring $\\big|1\\rangle$ in the last qubit, which corresponds to\n",
    "$\\sin^2(a*x+b)$.\n",
    "Together with the approximation above, this allows to approximate the values of interest.\n",
    "The smaller we choose $c_{approx}$, the better the approximation.\n",
    "However, since we are then estimating a property scaled by $c_{approx}$, the number of evaluation qubits $m$ needs to be adjusted accordingly.\n",
    "\n",
    "For more details on the approximation, we refer to:\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set the strike price (should be within the low and the high value of the uncertainty)\n",
    "strike_price_1 = 1.438\n",
    "strike_price_2 = 1.896\n",
    "strike_price_3 = 2*strike_price_2 - strike_price_1\n",
    "\n",
    "# set the approximation scaling for the payoff function\n",
    "c_approx = 0.25\n",
    "\n",
    "# setup piecewise linear objective fcuntion\n",
    "breakpoints = [uncertainty_model.low, strike_price_1, strike_price_2, strike_price_3]\n",
    "slopes = [0, -1, 1, 0]\n",
    "offsets = [1, 1, 1+strike_price_1 - strike_price_2, 1]\n",
    "f_min = 0\n",
    "f_max = strike_price_2 - strike_price_1\n",
    "butterfly_objective = PwlObjective(\n",
    "    uncertainty_model.num_target_qubits, \n",
    "    uncertainty_model.low, \n",
    "    uncertainty_model.high,\n",
    "    breakpoints,\n",
    "    slopes,\n",
    "    offsets,\n",
    "    f_min,\n",
    "    f_max,\n",
    "    c_approx\n",
    ")\n",
    "\n",
    "# construct circuit factory for payoff function\n",
    "butterfly = UnivariateProblem(\n",
    "    uncertainty_model,\n",
    "    butterfly_objective\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAE+CAYAAACazvcJAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd5xTdfb/8dcBRJoigiAowlrWvquLDcsqri62taEiYkVBLMjwdbu/Xdu6xS1DsQAKgopiwd5QEOzowq6uoigWQARRkCIMUobz++PcSAiZmWQmk09ucp6PRx5Dbu5N3hMyObd8iqgqzjnnXKYahA7gnHMuXrxwOOecy4oXDuecc1nxwuGccy4rXjicc85lxQuHc865rHjhcAVLRK4TEU26LRCRCSKyS8BMh4nIf0TkOxHRaFlzERkvIkuinBdWse2YlN8ncbszr7/Exjz9ROTUNMvniMg/QmRy8dAodADnarAcOC76987AjcBkEdlbVVcFyDMC+AroDqyJll0G/AI4H/gC+KSa7WcBF6Us+yrHGTPVD3gPeCxl+WnAkvzHcXHhhcMVuvWqOi369zQRmQe8ApwAPBQgzx7ASFV9KWXZh6o6IYPtVyX9PgVJVf8bOoMrbH6qysXNjOhnZwAR6SoiT4jIQhFZJSJvi0jvxMoism10WunC5CcR86mIlCctO1pE3ozWXyQit4lIi+ixo6JTUw2BIdEppjEiMge4GNg/ceqptr9Y4jVEZJ+U5VNF5OGk+2NEZLqIHCsi/4t+71dFZO+U7RqKyO9E5CMRWSMi80VkTOI5gS7ABUmnzC6MHtvsVJWInCUi70bP87mI3CQijZIevzB6jn1F5IUo0ywROb2274crXF44XNx0jn5+Gf3sBLyGfXn/ApgA3CUivQBU9RvgUeDClOc5CvgBMBog+tJ9DlgM9ACuBc4BEl/Y/wG6Rv/+Z/TvG7HTOs9gp6C6Jq1TJRFplHzL6Lfe3E7A34GbgF5AW+ABEZGkdUYA1wMPAicBVwPNoscujzI/k5T76Sry/hx4AHsPTgGGAb8Ebkmz+n3AE9j7MhsYLyI71vJ3dAXKT1W5gpf05bozcBvwLTAJQFXHJ60nwMvAjkBf4P7ooVHA8yKys6p+Gi27CJihqu9G9/8AzAVOVtXK6Pm+wb6Mu6rqG9ipMoA5yaebRORroF2Gp6C6AOtSfr/dVPXjDLZNti1wmKrOjp6jAVYgdwdmicgeWDEdqKpDk7Z7AEBV3xeRVcDXGeS+AZiqqhdE95+L3oe/iMifVHV+0rrlqpooxjOARVjRGp7l7+cKmB9xuELXGvuiXQd8iBWPnqq6EEBEWonIUBGZm7ReP+CHSc8xGSsKF0TbbIUdVdyVtM5BwKOJohGZAKwHDs/h7/MBcGDK7fNaPM+cRNGIvB/9TOzdd4t+jqnFc39PRBoCP2Hz60kPYN8fqUdYzyf+oapLsAv/fsRRZPyIwxW65cAxgGKnpxbopkM6jwEOwU4bvQ+swFo5nZJYQVVVRO4C+ojIdcBZ2LWK+5Kepz22d0zSdpUisgTbu8+VClWdnoPnWZZyf230s0n0szV2IX5FHV+nDbAFKe9N0v3U9yZdria4ouKFwxW69VV90YpIE+w0yBWqOjxpeboj6buw6xbdsOsdj6nq0qTHF2LXCZKfvyH2BfxNXX6BLHwX/WycsrwVdu0lG0uA5iKydR2Lx2LsKK5tyvJ20c98vTeugPipKhdnW2Kf4UR/isRpqJNTV1TVz7HTKNdjp57uSlnlTeC0qFgknI7tXL2a29hVSlwr2DOxQEQ6Ys19s/Vi9PP8atap8WggOnU3Azgz5aGzgA3AG7XI5mLOjzhcbKnqchH5N/BHEVmBfZH9Fju9tXWaTUZh5+rnAy+kPPYn4L/AYyJyO3Ze/m/AxOjCeL1T1fkiMh24UUQqsKL4e2qxV6+qH4rISOCfItIWazSwDXCGqp4drTYL6C4i3bEjlM+i6xKprgUmRqf7xgP7YqcG70i5MO5KhB9xuLg7B/gUuBsYgl3QvruKdZ/CLnaPVdUNyQ+o6kzgeOyUzCNYIbkfOKN+YlepFzAPuBf4M9ai6cNaPtfl2BHWuViz28FARdLjf8Iu1j8I/BtrzrwZVX0eOBs4AHgSKMOaJF9Zy1wu5sSnjnWlQkROwIrHD2vR/NU5F/HC4YqeiHQAdsM6rs1T1ZMCR3Iu1vxUlSsF/bC+HN8BAwJncS72/IjDOedcVvyIwznnXFZKojlumzZttHPnzrXadtWqVTRv3jy3gepRnPLGKSvEK2+cskK88sYpK9Qt74wZMxar6nabPaCqRX/r0qWL1taUKVNqvW0Iccobp6yq8cobp6yq8cobp6yqdcsLTNc036l+qso551xWvHA455zLihcO55xzWcl74RCRXUVkRDTlZWU0hWUm27UUkbtEZKmILBeRcSLSup7jOuecSxGiVdXewAnANGyc/0w9iE3Ocwk2mN3fgMeAI3Id0DnnXNVCnKp6UlU7quqZwMxMNhCRrsDPgQtUdYKqPooN3Ha4iBxTLynHjYPOnTny6KOhc2e7X8jiltc5iNfnNk5ZoV7z5v2IQ1NGJc3Q8cAiVX056XneEpHPoscm5SofYG9wv35QUYEAzJ1r9wF6987pS+VE3PI6B/H63MYpK9R73rh0ANwDmzsg1QfUbpKb6l1zDVRUbLqsogKuuAI+rO0I1/Vo6ND0ea+5pjA/1M5BvP7OqvobK8SsUO/fCUHHqhKRh4E2qnpUDeu9gM2ffGrK8nuBnVX10DTb9MMGt6Ndu3Zdxo8fn3GuI48+Gknzvqg9ccbPkzeqpEulIrz04otpHikMK1eupEWLFqFjZCxOeeOQNVZ/Z1X9jUHhZYWcfSd069ZthqoekLo8LkccWVPVkcBIgAMOOECPOuqozDfeaSc7tEshnTrBnDm5CZhLnTunz7vTTmT1e+fZ1KlTCzpfqjjljUXWjh1h3rzNFhfk31lVf2OFmBXq/TshLv04lgIt0yxvFT2WWzfdBM2abbqsWTNbXojiltc5gPPO23xZoX5u4/Y3Vs9541I4ZpH+WkZV1z7qpndvGDkSOnVCRaBTJ7tfqNcLkvMCNGwIt99euHmdA3jvPWjRAjp2LPy/szh/J9RD3rgUjmeB7UXk8MQCETkA2Dl6LPd694Y5c+x84Jw5hfsBSYjyvnfTTVBZCU2ahE7kXNU+/hieeALKymDevHj8ncX0O6E+8oboOd5MRM4QkTOAHYDtEvdFpFm0zsciMiqxjaq+ATwP3C0ip4vIqcA44FVVzW1T3JhbcsghsOuuUF4eOopzVRs6FBo1gssvD53E1UKII462wEPR7RBgr6T7baN1GgENU7brCbwEjAbuBmYAp+Uhb7w0aAADB8K0aXZzrtAsWwajR0OvXtC+feg0rhbyXjhUdY6qShW3OdE6nVX1wpTtlqnqRaq6japurarnqOrifOePhQsvhJYtYfDg0Emc29ydd8KqVXaaysVSXK5xuGy0aAF9+8LDD8Pnn4dO49xG69fDsGFw5JGw//6h07ha8sJRrK68ElThlltCJ3Fuo0cftb4bgwaFTuLqwAtHserUCXr0sCZ4q1aFTuOcKS+HXXaBk04KncTVgReOYlZWZhcix44NncQ5ePNNeOMNa7zRMLXti4sTLxzFrGtXOOggGDIENtRmUGLncqi8HLbe2hpvuFjzwlHMROyo46OP4Nn66SfpXEbmzbPGGn37wlZbhU7j6sgLR7E74wzYYQdvmuvCuuUWa6wxYEDoJC4HvHAUuy22sBZWkybBu++GTuNK0cqVcMcd1lijU6fQaVwOeOEoBf36QdOmdq3DuXwbO9YaaXgT3KLhhaMUbLstXHAB3HsvfPVV6DSulGzYYDssBx9sjTVcUfDCUSoGDoQ1a2DEiNBJXCl5+mmYPduHFykyXjhKxR57wPHHw623WgFxLh/Ky2HHHe36hisaXjhKSVkZLFoEDzwQOokrBW+/DVOmWEuqLbYIncblkBeOUnLssbDXXrYXqBo6jSt2Q4bYdKV9+4ZO4nLMC0cpSXQIfPttePnl0GlcMfvyS7jvPrjoImjVKnQal2NeOErNuedC69beIdDVr9tvh7VrrVGGKzpeOEpN06bQvz88/jh88knoNK4YffedFY6TToLddgudxtUDLxyl6PLLbb7nYcNCJ3HFaNw4+Ppr7/BXxLxwlKIOHaBnTxg1CpYvD53GFRNVa3zxox9Bt26h07h64oWjVJWV2RhCo0eHTuKKyeTJMHOmHW2IhE7j6okXjlLVpQsccQQMHQqVlaHTuGJRXg7t2kGvXqGTuHrkhaOUDRoEc+bYhXLn6mrWLHjmGbuGtuWWodO4euSFo5SdfDJ07mx7ic7V1ZAhVjD69w+dxNUzLxylrGFDuOoqePVVmD49dBoXZ0uW2PDpvXtD27ah07h65oWj1F18sU3l6R0CXV2MHAmrV/souCXCC0ep23pr6NPHBj5csCB0GhdHa9fa1LDHHAP77hs6jcsDLxzOTldVVtqQ685l6+GHbafDO/yVDC8cDnbeGU45xSZ5qqgIncbFSaLD3+67w3HHhU7j8sQLhzODBtkFznvvDZ3Exclrr1nDioEDoYF/nZQK/5925ogjYP/97SK5z9XhMlVebsOmn39+6CQuj/JeOERkLxGZLCIVIrJARG4QkYYZbLe3iDwfbbdYRG4XkRb5yFwSROyo44MP4PnnQ6dxcfDZZ/DYY3DppdC8eeg0Lo/yWjhEpBUwCVDgFOAG4Grg+hq2awm8CDQFegK/BHoAfl4ll3r2hO2396a5LjNDh9rpqSuvDJ3E5VmjPL9ef+zL/3RVXQG8ICJbA9eJyM3RsnQuj7b7haouAxCRJcATInKAqnrvtVxo3BiuuAL+8Ac78thzz9CJXKFascJGVz7rLNhhh9BpXJ7l+1TV8cDElAIxHisKR1az3X7A9ETRiLyAHbmcmPOUpezSS23YiCFDQidxhWz0aPj2W+/wV6LyXTj2AGYlL1DVeUBF9FhVmgBrU5atBzYAvlucS9ttB+edB3ffba2snEtVWWmnqQ47DA48MHQaF0C+T1W1ApalWb40eqwqHwPniMgWqrouWtYFaAhsm24DEekH9ANo164dU6dOrVXglStX1nrbEHKRt/mhh3LgnXfy6W9/y7zevXMTLI1SfG/zpT6ztnn5Zfb57DPeu+ACFufoNfy9rT/1kldV83YD1gFlaZbPB/5czXZ7YEcYw4Htgb2BGdGy52p63S5dumhtTZkypdbbhpCzvMceq9qhg+ratbl5vjRK9r3Ng3rNevjhqp07q65fn7On9Pe2/tQlL3aJYLPv1HyfqloKtEyzvFX0WFqqOgs7eugFLAT+B7wFvA18mfuYjrIyG0bioYdCJ3GFZPp0G035qqtsdGVXkvJdOGaRci1DRDoCzUi59pFKVUcD7YAfAR2AK4FdgWn1krTUHXecDSNRXu4dAt1GgwfbaMoXXxw6iQso34XjWaC7iGyVtKwnsBp4qaaNVfU7VX1XVRcB52L5H6yXpKWuQQMbRmL6dHj99dBpXCH44gsbRblPHxtV2ZWsfBeO4cAa4BEROSa6gH0d8C9NaqIrIh+LyKik+1uLyN9E5EQR6S4ifwXuBK5S1W/y/DuUjvPPt+EkvEOgAxs9ubLSTlO5kpbXwqGqS4GfYa2hnsR6jJcD16as2ihaJ6ES2B+4B3gseo4zVXVMPUcubc2bQ79+8MgjNje5K10VFTZ68qmn2mjKrqTlfawqVX1fVY9W1aaq2l5V/6CqlSnrdFbVC5Pur1LVn6vqttF2B6rqY/nOXpKuuMLGsbrlltBJXEh33w3ffONzbjjAR8d1NenYEc48E+64w3oKu9KzYYOdruzSBQ4/PHQaVwC8cLialZXZ2ERjxoRO4kKYOBE+/NCONkRCp3EFwAuHq9nBB0PXrjZ+VWVlzeu74lJeDu3b25Gnc3jhcJkqK4NPPoGnnw6dxOXTe+/BCy/Y0OmNG4dO4wqEFw6XmdNPt+sd5eWhk7h8GjwYmja1UZOdi3jhcJlp1AgGDICpU+Htt0Oncfnw1Vc2B/3550Pr1qHTuALihcNl7pJLrG+HdwgsDcOHw5o1PueG24wXDpe5Vq3gwgvh/vvhSx9bsqitWQO33QbHHw97VDdVjitFXjhcdgYOhLVr4fbbQydx9Wn8eFi0yI82XFpeOFx2dtsNTjrJCsd334VO4+qDqjWC2HtvOPbY0GlcAfLC4bI3aBB8/bWdsnLFZ+pUeOcdO9rwDn8uDS8cLnvdusG++/pcHcWqvBzatIF6nDbYxZsXDpc9EdsbffddmDIldBqXS7Nnw1NPwWWXWf8N59LwwuFq55xzYLvtvENgsRk61PrsXH556CSugHnhcLXTpIntlT71lO2luvhbtgzuugt69YLttw+dxhUwLxyu9i67zMYvGjIkdBKXC3fcAatW+ZwbrkZeOFztbb+97Z3edRcsXRo6jauL9eth2DA46ijYb7/QaVyB88Lh6qaszKYVvfPO0ElcXUyYAJ9/7kcbLiNeOFzd7Lef7aUOG2Z7rS6eysthl12sc6dzNfDC4epu0CDbW33kkdBJXG1MmwZvvmnDyTTwrwRXM/+UuLo78UTbW/VRc+OpvBxatoSLLgqdxMWEFw5Xdw0b2t7qG2/YnquLj3nz7PpG377QokXoNC4mvHC43LjwQth6az/qiJthw+zngAFhc7hY8cLhcmOrrWyv9aGH7HqHK3wrV1rfjR49YKedQqdxMeKFw+XOlVfaoIe33ho6icvEXXfB8uXeBNdlzQuHy53OneH002HkSOuB7ArXhg3W4//gg+GQQ0KncTHjhcPlVlmZ9SK/++7QSVx1nnoKPvnEjzZcrdRYOETkfBFpnY8wrggceigceKBdJN+wIXQaV5XycujY0a5vOJelTI447gJ2ARCRShE5qH4juVhLzNXx0Ufw3HOh07h03n7bZvkbMMCGUHcuS5kUjqVAh+jfAviUb656Z5wBHTr4XB2Fqrwcmje3VnDO1UImhWMScI+ITMeKxhgReauqW01PJiJ7ichkEakQkQUicoOINMxguwNE5HkR+Sa6TRKRgzPI7/KtcWNrYTVpErz3Xug0LtnChTZX/EUXwTbbhE7jYiqTwtEHuB74L3bE8Rkws5pblUSkFVaIFDgFuAG4Onr+6rbrGG3XCDgvujUCXhCRThn8Di7f+vWzqUd9ro7CcvvtNhjlVVeFTuJirMYTnKpaAfwDQESOAa5R1Xdq+Xr9gabA6aq6Avvi3xq4TkRujpalcyKwFXCaqi6PsrwOLAZOAG6vZR5XX1q3hvPPhzFj4M9/tmlmXVirV1vhOOkk2G230GlcjGXSqqpSRA6M7k4Fqvpyz8TxwMSUAjEeKyZHVrPdFsB6ILlzwMpomdQhj6tPAwfCmjUwYkToJA5g3DhYvNib4Lo6y+RU1Vpgy+jf5wN12XXcA5iVvEBV5wEV0WNVmRCt808RaSsibYFy7ML9Q3XI4+rTnnvCccdZT/I1a0KnKW2q1kT6xz+2+VOcq4NM2uK9j51Kegzbuz9DRA6oYl1V1epOG7UClqVZvjR6rKonXSAi3YCngMTJ2YVAd1X9uqZfwAVUVmbF48EH4bzzQqcpXS+8ADNn2qlD8YN0VzeiWn3rWhE5FBiBHRE0oPpTQ6qqVbaQEpF1wK9UdXDK8vnA3ar6+yq2aw+8jBWxRGG6AtgfODQ6akndph/QD6Bdu3Zdxo8fX03sqq1cuZIWMRpuuuDyqnLgRRexoXFjZowYscmXVsFlrUGc8qZm3fc3v2Gr2bN5Y/x4tHHjgMnSi/N7W+jqkrdbt24zVHXzAwVVzfgGbAAOymablO2/Aq5Ns3wVVlCq2u5fwBxgi6RljYG5wNCaXrdLly5aW1OmTKn1tiEUZN4RI1RB9aWXNllckFmrEae8m2R9/317/6+/PliemsT2vY2BuuQFpmua79Rsx6rqhu3119YsUq5lRE1tm5Fy7SPFHsBMVV2XWKCqa7Hmv7vUIY/Lh/POs1ZW3iEwjCFDYMstoX//0ElckciqcKjqS6q6UkQOFpGrReSm6GemHfGeBbqLyFZJy3oCq4GXqtluLrCPiHx/jC0iWwL7YEcirpA1bQqXXgqPPw6ffho6TWlZssQGnDz3XGjbNnQaVySyKhwi0lxEngFeB/6CdQ78C/C6iDwtIs1qeIrhwBrgERE5JroOcR3wL01qoisiH4vIqKTt7sSGPXlURE4UkZOAx4D2wMhsfgcXyBVX2BSzQ4eGTlJaRoyw/htlZaGTuCKS7amqm4GuwNlAE1VtDzSJ7ncF/lbdxqq6FPgZ0BB4EusxXg5cm7Jqo2idxHYzgOOwToD3AHdjp7eO1dp3RnT51KED9OwJo0fDirp0BXIZW7sWbrkFjj0W9tkndBpXRLIdGrMH8BtV/b7vhKpuAB6KhhO5Aah28mJVfR84uoZ1OqdZNhmYnGVeV0gGDbJOaKNGeSe0fHjwQRub6s47QydxRSbbI46WQFUTSn8ObF23OK6odekChx9up6sqK0OnKW6q1hhh992tH41zOZRt4XgHuExk0x5E0f3Losedq9qgQTBnjl0od/Wm5bvvwn/+Y9c2GvhEny63sj1V9XusZdQsEXkUWAS0BU4DOmNjUTlXtVNOgTZt4JxzOHLtWthpJ7jpJujdO3Sy4jBuHFxzDfvNnWsFY8sta97GuSxlVThU9UUR+QnwB+BMrFXTQuBNbMTbuvTxcKVg/HhYvhzWrbMhCObOtSHYwYtHXY0bZ+9lRYW9txs22LwojRv7e+tyKutjWFWdqapnq+ouqtos+nmOFw2XkWuugXXrNl1WUWHLXd1cc429l8n8vXX1INt+HP8Ukb3qK4wrAfM2G1as+uUuc/7eujzJ9ojjNODdaJrY/iLSsj5CuSK2007ZLXeZ8/fW5Um2Q47sDByDjSv1D2ChiNwXzQzoXM1uugmapQww0LSpLXd1c+ONmw+Z3qyZv7cu52pzjWOKqp4PbI919tsRmCgic0XkehHZOdchXRHp3RtGjoROnfh+QP9evfzibS60aGH9N9q0QUWgUyd7r/29dTlW6wbeqrpSVUdhw4W8BnQEfgd8JCKPi0inHGV0xaZ3b5gzh5defBH23x/eeMO+8FzdlJdbsVi40N7bOXO8aLh6UavCISKdReRaEfkUeB6b//tMbCypk7E+HbWbOcmVDhHrEPjBBzBxYug08TZjBrzyClx1FTTKtnuWc9nJtlXV+SLyIvAxcAFwF/ADVT1BVSeo6hpVfQab3rWq6WWd26hnT2jf3ubDdrU3ZIidqrr44tBJXAnI9ohjBPAlNtf3zqp6o6rOT7PeR8Cf6pzOFb/GjW3I9YkT4X3vClQrCxdax8qLLoKW3tDR1b9sC0eHqLNftaPUqupCVb2+DrlcKbn0UmjSxI86auu222D9ejtN5VweZNscd2l9BXElrE0bm172nntg8eLQaeJl9WoYPhx+8QvYddfQaVyJyPriuIj0FJFJIjJPRL5KvdVHSFcCysrgu+9sxjqXuXHjrNj6/CYuj7K9OH4OMBa7OL4j8ATwVPQ8K4Bbch3QlYi99oLu3eHWW23mOlczVTu99+Mfw5FHhk7jSki2Rxy/Am4Eroju36aqfYAfAIuBiqo2dK5GZWV2ofeBB0IniYdJk2DmTHvfUnuMO1ePsi0cuwGvqWolUEk045+qfovNN35lbuO5ktK9O+y5p3Vk8w6BNRs8GNq2tZ73zuVRtoVjBZCYGeYLYM+kxwRonYtQrkSJ2N7zf/8LL78cOk1hmzULnnkGLr/cJ2tyeZdt4fg38KPo308AfxSRviJyAfB3YFouw7kSdN550Lq1N82tydCh1gfmsstCJ3ElKNvC8RcgMbj/H4G3gNuxHuSLgUtzF82VpKZNrV/H44/DJ5+ETlOYvvkGxo61cajatg2dxpWgjAqHiDQVkR7AYUAjEWmnqstU9RSgObCNqh6sqp/WZ1hXIq64wsZbGjo0dJLCdMcdNrPfwIGhk7gSVWPhiIZJnwk8hJ2Ougf4UER+DhCNT7WiXlO60tKhg41hNXq0zU/uNlq3Dm65BY4+2prhOhdAJkccNwMbgCOAZsDewH+xcaucqx+DBsHKlTBqVOgkhWXCBJg/3xoROBdIJoWjK/D/VPU1Vf1OVT/ArmXsJCLt6zeeK1k/+Qn89Kd2umr9+tBpCsfgwTa0yIknhk7iSlgmhaM9kHrt4hOs+e32OU/kXEJZGcydC489FjpJYZg2Dd58065tNKj1HGzO1Vmmnz7vjeXy7+STYeedrUOgs/ehZUu48MLQSVyJy7RwTEwZyHBhtHyyD3Lo6k3DhjZU+Ouvw1tvhU4T1rx5dn2jb1+bsMm5gDKZY9Ln1XDh9OkDf/yjndu/777QacK55RYbhuVKH9XHhVdj4fAJmVxQW21l06EOGwY33ww77hg6Uf6tXGl9N3r0gE6dQqdxLvv5OOpKRPYSkckiUiEiC0TkBhFpWMM214mIVnH7Xb6yu0Cuugo2bLC97lJ0992wbJk3wXUFI6+FQ0RaAZOwi+2nADcAV1Pz6bA7sWbBybe/RY89Wy9hXeHo3BlOOw1GjoRVq0Knya8NG+w03UEHQdeuodM4B+T/iKM/0BQ4XVVfUNXhWNH4PxHZuqqNVHW+qk5LvgH7ArNU9e38RHdBDRoES5fa3ncpefZZmD3b59xwBSXfheN4YGLKECXjsWKS8RRmItIaOBa4P7fxXME69FA44ADb+96wIXSa/Ckvhx12gDPOCJ3Eue/lu3DsAcxKXqCq87CZA/fI4nl6AFvghaN0iNhRx0cf2V54KXj3XZg82VpSbbFF6DTOfU80jzOticg64FeqOjhl+XzgblX9fYbP8yLQUlW7VLNOP6AfQLt27bqMHz++VplXrlxJixi1m49T3myzyvr1HNKrFxU77cQ7//xnPSZLL9/v7e5//zttJ0/mjQcfZP3WVZ7JTStOnwOIV944ZYW65e3WrdsMVT1gswdUNW83YB1Qlmb5fODPGT5He2za2l9m+rpdunTR2poyZUqttw0hTnlrlfUvf1EF1f/9L+d5apLX945Z5TQAAB+BSURBVHbRItUtt1Tt379Wm8fpc6Aar7xxyqpat7zAdE3znZrvU1VLgZZplreKHsvEWdg4WQ/kKpSLkX79bLKnYp8hcMQIWLPGmiI7V2DyXThmkXItQ0Q6YsO1z0q7xebOBl5V1c9znM3FwbbbwgUXwLhx8FWRjnCzZg3ceiscdxzsuWfoNM5tJt+F41mgu4hslbSsJ7AaeKmmjUWkM3AIflG8tJWV2Zfr7beHTlI/HngAFi2yxgDOFaB8F47hwBrgERE5JrqAfR3wL01qoisiH4tIuhl8zgbWY7MRulK1++5wwglw221WQIqJqp2G22svOPbY0GmcSyuvhUNVlwI/AxoCT2Kd/8qBa1NWbRStk+psYLKqLq7PnC4GysrsVNX9RXbw+fLL8N//eoc/V9AyGR03p1T1feDoGtbpXMXy/eojk4uhY46BffaxDnIXXFA8X7KDB0Pr1nDuuaGTOFcln0bMxZOI7ZX/738wZUroNLnxySfw+ONw6aXWcsy5AuWFw8VX796w3XbF0zR32DCbvOqKK0Inca5aXjhcfDVpApddBk89ZQMBxtmKFTB6NPTsCR06hE7jXLW8cLh4u+wyG8dpyJDQSepm1Cj49ltvgutiwQuHi7ftt4deveCuu2zY9TiqrIShQ+Hww6FLlcOvOVcwvHC4+Bs0CCoqbHrVOHr8cZgzx2f4c7HhhcPF349/DN262cXldetCp8ne4ME2y+Gpp4ZO4lxGvHC44lBWBvPnwyOPhE6SnRkz4JVXYMAAa1HlXAx44XDF4aSTYNddrUNgnAweDC1awMUXh07iXMa8cLji0KABDBwIb74Jb7wROk1mFiywAQ379IGW6WYbcK4weeFwxePCC+0LOC4dAm+7Ddav9zk3XOx44XDFo0UL6NsXJkyAefNCp6ne6tUwfDicfDLsskvoNM5lxQuHKy4DBtjPYcPC5qjJuHGwZIk3wXWx5IXDFZeddoIePaxPx8qVodOkl5hzY7/94MgjQ6dxLmteOFzxGTQIli+HMWNCJ0lv0iSYOdPn3HCx5YXDFZ9DDoGDD7bxqzZsCJ1mc+Xl0K4dnH126CTO1YoXDlecBg2Cjz+2kXMLyaxZ8OyzcPnlsOWWodM4VyteOFxx6tEDOnYsvA6BQ4dawejfP3QS52rNC4crTo0aWQurqVPh7bdDpzHffANjx9oEVG3bhk7jXK154XDF65JLoHnzwukQeMcdNoqvN8F1MeeFwxWvVq2sN/n998OXX4bNsm6d9S352c9g333DZnGujrxwuOI2cKB9ad92W9gcEybAF1/40YYrCl44XHHbbTcbOff2222Yj1AGD7YsJ5wQLoNzOeKFwxW/QYNg8WK4774wr//GGzZq78CBNoqvczHnn2JX/I46Cn70I2uaq5r/1x88GLbZBi64IP+v7Vw98MLhip+IHXXMnGnDfeTTvHl2faNvXxu917ki4IXDlYZevWyYj3x3CLzlFvt55ZX5fV3n6pEXDlcattzShvl49lkb9iMfVq60vhunn26j9jpXJLxwuNLRv78VkCFD8vN6Y8fCsmV2msy5IuKFw5WOtm1tuI+xY20Spfq0YYMVqIMOstF6nSsieS8cIrKXiEwWkQoRWSAiN4hIwwy3PV1E/i0iq0VkiYg8JyLN6zuzKyJlZdafY+TI+n2dZ56B2bPtaMPn3HBFJq+FQ0RaAZMABU4BbgCuBq7PYNtLgPuAZ4HjgUuA2UCj+srritC++8Ixx9hF63Xr6u91Bg+GHXawUXqdKzL5PuLoDzQFTlfVF1R1OFY0/k9Etq5qIxFpA5QDA1T1j6o6VVUfVdUBqro8P9Fd0SgrgwUL4KGH6uf5330XJk+2llRbbFE/r+FcQPkuHMcDE1V1RdKy8VgxqW7y5bOin2PrK5grIccfD7vvXn8dAgcPhqZNoV+/3D+3cwUg34VjD2CTtpCqOg+oiB6rysHAh8DFIjJfRNaJyJsicmj9RXVFq0EDG/5j+nR47bXcPvdXX8G4cdZLfNttc/vczhUI0TwOwSAi64BfqerglOXzgbtV9fdVbDcROBRYAfwaWBL9PADYTVUXpdmmH9APoF27dl3Gjx9fq8wrV66kRYx6/MYpb8isDVavpmvPnizbf39mXl/jJTYgs7ydxo7lB2PG8NbYsVQE7LsRp88BxCtvnLJC3fJ269ZthqoesNkDqpq3G7AOKEuzfD7w52q2ex67oH5c0rKtgaXAjTW9bpcuXbS2pkyZUuttQ4hT3uBZf/tb1QYNVD/9NKPVa8z73Xeq7dqpHn983bPVUfD3NktxyhunrKp1ywtM1zTfqfk+VbUUaJlmeavoseq2U2BqYoHadZIZwF45zOdKyRVX2GmrYcNy83wPPACLFnmHP1f08l04ZpFyLUNEOgLNSLn2keIDQKLbJpsDG3IZ0JWQHXeEM8+EO++EFStqXr86qnaxfa+9rLmvc0Us34XjWaC7iGyVtKwnsBp4qZrtnop+dkssEJGWQBfgnVyHdCVk0CD49lsYPbpuz/Pyy/D229bU1zv8uSKX78IxHFgDPCIix0QXsK8D/qVJTXRF5GMRGZW4r6rTgceBUSJygYicCDyBXTO5NZ+/gCsyBx4Ihx0GQ4dCZWXtn6e8HFq3hnPPzV025wpUXguHqi4FfgY0BJ7EOv+VA9emrNooWifZucBjwL+Ah7GicXT0nM7VXlkZfPYZPPFE7bb/5BPbtn9/67/hXJHL+3Adqvo+cHQN63ROs2wlcFl0cy53Tj0VOnWyo4bTTst++2HDoFEjG7bduRLgo+M616gRXHUVvPIKzJiR3bbLl8OoUdCzJ3ToUD/5nCswXjicA7j4YpvadfDgmtdNNnq0TdhUVlY/uZwrQF44nANo2RL69IHx420AxExUVtpF9SOOgC5d6jefcwXEC4dzCVddZcXg1gwb6j3+OMyZ40cbruR44XAuYZdd4JRTYPhwqKioef3Bg6FzZ9vGuRLihcO5ZIMGwTffwL33Vr/ejBl2Mf2qq6BhRhNYOlc0vHA4l+yII2D//e1oorqRowcPtovpffrkL5tzBcILh3PJROyo44MPYOLE9OssWGAX0S++2C6qO1divHA4l6pnT2jf3joEpnPbbXYRfcCA/OZyrkB44XAuVePGNuT688/DzJmbPrZ6tV08P+UUu5juXAnywuFcOpdeCk2awJAhmy6/915YssSb4LqS5oXDuXTatIHzzoN77oHFi22Zql0U328/+OlPw+ZzLiAvHM5VpawMvvvOTk0BraZPh/fft4vnPueGK2FeOJyryl57Qffu1pN8zRp2nDAB2rWzi+fOlTAvHM5VZ9Ag+PJLaNuW1m++aUcgDz8cOpVzQXnhcK46ixfbaanEnOTLl0O/fjBuXNhczgXkhcO56lxzzeY9yCsqbLlzJcoLh3PVmTcvu+XOlQAvHM5VZ6edslvuXAnwwuFcdW66CZo123RZs2a23LkS5YXDuer07g0jR0KnTqgIdOpk93v3Dp3MuWC8cDhXk969Yc4cXnrxRZvxz4uGK3FeOJxzzmXFC4dzzrmseOFwzjmXFS8czjnnsuKFwznnXFZEU4dTKEIi8jUwt5abtwEW5zBOfYtT3jhlhXjljVNWiFfeOGWFuuXtpKrbpS4sicJRFyIyXVUPCJ0jU3HKG6esEK+8ccoK8cobp6xQP3n9VJVzzrmseOFwzjmXFS8cNRsZOkCW4pQ3TlkhXnnjlBXilTdOWaEe8vo1Duecc1nxIw7nnHNZ8cLhnHMuK144nHPOZcULh3POuax44XDOOZeVRqEDuNoRkY7ACYAAD6nqEhHZEfglsAswBxipqu+GSwki8hvgmdA5MiUiTYFGqvpt0rLtgCuBvYANwNvAbaq6PExK58Ly5rgRERHgNOBEYE9gW+xL4ktgGjBGVT8Kl3AjETkIeAFoDqwHvgG6A88AlcBMYB9ge+AYVX0lUFREZAOgwCzgPmC8qn4SKk9NROQZYLaqDozudwWexT4LM7BC3QVYCxytqjMDZt0faKqqryctOw74HRuL3DvAdcnrFIrob+4XwE+wz8h0bCejoL+URGRrbOyno1X11dB54PtMRwONgadVdVW0w3MFsDPwKbYjuSAnr1fg/0d5Eb3Bz2BfCF9iXwo7YB/mZ7E3fnfgRlW9MVTOBBF5ATtaPA1YBZQDp2JfbGeo6joR2RJ4DGiiqt0CZt0A/A3YFzgWy/0frIg8qKpfhMqWjogsBi5W1cej+9Ow9/jUxFGIiLQEngC+U9XuAbNOA55U1Zui+32AO4EpwItYkfsZcATQI/E7Bcr6Ova+fhDdbwU8j/3NrYxWa4HtpHVPPuILQUQur+bhpsDfgSHAbABVvS0fudIRkV2ByUDHaNFnwM+xncttgE+w76/VQBdVnV/nF1XVkr8B92MfgH2TlnUAngMmRPePxD7gfQog7xLg+KT7bbG9y5+nrHcisDhw1g3AQdG/twX6RR/y9dFtarRs29Dva5SxAvhp0v21qe9r0nu7KnDWFcnZgI+BYWnWGw68Uyifg+j+KOxI+bikZccBS4HyAvgcbMCO3jdUcUt+rDJw1gexI8tdo7+xe6Lvs9eBraJ12kTrjMjFa/rFcXM88FtNOg+vdkjXHzhVRNqr6kvAn4GBgTKm0jT/Tj18LKjDSVX9RlVHqurPgB2Bq7FD6+HAQhF5OmhA8x6QfIS2CPtjTNUaKzIhbUi53wl4OM16D2N7nIXkZOAGVX0usSD6903A6cFSbfQ48BVwMdBQVRskbtjnQYCjomUNQwYFDgduUtWPVfUb4P9h1zn/odGRm6ouBgaz6We71rxwmAbYHkSqSuwD0jK6/ybww3yFqsZ04Jci0kJEGgC/B74ALhORhgAi0gi4HPsiLDiq+qWqDlHVQ4EfANdiR3mh/RX4rYj0id7Dm4C/i8ixItJYRLaMriP8BdvTC+kVoHfS/ZlAuuGzD8Q+H4VkG+xznGoGdm0uKFU9DbgA+BXwbxE5LPnhMKmq1Ao7xZ6Q+L9OnYPoU2yHrc68VZV5AfiTiPxPVT+F78/BDsX+QxIXxVsAhdCS5hos81LsdE8FdmHsIWC2iCQujnfADv8LmqrOxb6w/1oAWR4RkQHY3lk58CG24/Ac9oUh0apPYF8qIf0eeC3aeRiGXRQfKyLbYqcABftclAG/DRUySQ8RSRS2pcBmEwRhp1RW5C9S1VT1eRH5Efb+PS0iz2GtFoNef0njK+xoM6ESGIEdLSdrS46y+8VxIGrGOhE7mpiLndf+AbAG6KWqz0br3YzNiNUzVNaEKPNJWPGfoKoLRWR74NfYaYm5wJ2q+p+AMRGRa4E7NEetOfJFRFoDPYGDsD1gwb7sPgCeUtUZAeN9T0T2A24HDmbTwpb491LslNCQMAlN1Egi1RhV7ZOy3ghgL1U9Ij/JMhP9bd2MnUYbgRWTbqr6ctBggIg8BnyT+l6mWW8YsIeqHlvn1/TCYaJTPGcBPwaaYBca74vOGTpX0ERkT6x4pBa511V1Xchs2RCRvsAnqvpi6CzpRM2zy7GdsxO1AJo5i0g7oJmqflbDev+HNZKYXOfX9MJRfESkgaqm28MrGCLSBLuAtwH4uBC/3KJrHDuT1KdHVeeFTeVceH5xPIWI7C0iPUTkkujWQ0T2Dp0rlYicLiKPicgzIvKLaFlPEZkDrBORudHeW1Aicm7UvyBxv5GI/BVrivk/7OL9NyJSCOfgARCRLiLyBHY++APgNeAN4DMR+UJEbhCRZkFDFhGJhM6Rjog0Tf2/FpH9ou+FLqFyBRey/XEh3YA+2HWBdG23K7EhPC4KnTPKelaU61Ws2WAF0Be7NjMK6y16f5S7e+Cs7wOXJd3/Z5T3D8BhWFPC67DOSb8vgPf259i1relY8+vrsE6ga6PMV2Otl94GWhVA3pOwfjEfRJ+Fn6ZZ52DC9zX4OVGfgqRlp2KdQdcD66L3/MTQ72mUrSXwaJRrPXAH0BAYm/K98CrQJnTeDH+nHrn6HAT/ZQrhBgyIPiC3Yr1s20QfkobRvw8Hbom+UK4ogLz/BoYn3e8dZftnynp3AZMCZ60Ajky6/xUwMM16vwTmFsB7OwMYW8VnZA52lN4k+sK7LXDWY6Mvr9eiz+eM6P4/iU5DR+sVQuGoZNMOgKdFX76vR//3V2NHdetJ0+EyQN6h2LAiA4Dzo52FCcDnURHcDuv/9QVwe+i8Gf5OOSscfo0DEJFPsS/im2tY79dAf1XdOT/JqsyxAjhdVSdF91tiF0OP0aSLitEprBGqGqx/hIgsBK5U1QnR/TXYUdDUlPWOBZ5Q1ab5T7lJjtXAyar6QsryVliP/b1V9QMROR/4m6q2D5EzyvQqNq7WRUnL+mBfei9gLQK/E5GDsYvkwTqqRa2qDlHVt6L7/wG+UNVfpKz3DNBcVY8MEDM5xxysU90d0f39scJ8kaqOTVqvL3ak/IMgQS3D6AxX7YR1Wqzz58CvcZj2wFsZrPcWBdA5CWtqmfyfnxjrZ1nKeiuxjlYhPYF1Vmwc3Z8E9EqzXi9sry60r7CWdal+jL3viX48c9nYMTSUfYB7kxeo6mhseJxDgBejPh2FaB+sWWuqkdigh6Ftx8b+WxCNSYWN+5TsY9L3R8mnC7CjoH1ruHWq6gmy5R0AzTtAXxF5WatojRRdvOuLXdANbS42+ulEAFWtjJoJfpCy3i7AwjxnS/U7rIfzeyJyJ/Ak8DcR2YeNndS6AftjI6WGNhK4UURaYEVuLdbz+hpgim7sj7IzELqF1XfYCMmbUNUZUU/nidipoOvynKsqyac3lpO+M9oqCmOH9jOsAL8U3T8CO7V2KHZdI+Ewwn8OZgNvqer51a0kImcAD+TiBb1wmKuxnsHvi8gj2BDgib33lsAe2DnZHSmMntiPkDJ0gKq+mWa9s9n0Q553qvqNiByCffH+H9Z7FaBrdFuLnVY5QlX/HSblRqp6U3Ra5bfAHxOLscYGZUmrrsMunof0P+w8+xOpD6jqp1HxeAYYk+dcVZkoIuujf7fEdhZeSllnD8Lv7ICNnzZERPbFitxZ2E7QH6KdinewI6NBQOgRs6dhBa0myR1E68SvcUREZBes1/VxbByeOOFzrGXN37WA55JIJSI7ActUtSCGcAAQkc5s2kntEy3MPhxbYEdsTYBPC+k9TBCRS7FhR/bXKjqqikhzrHXQMWoD9AURjSCQaraq3pey3tRoeSE0Jb8KO4W6BTYKw3AR6YVdQ0oMcjkS+E3Iz3DULPgwVR1aw3ptsGt0qcU6+9f0wrG5qN124trAMlUNPQqqc65ARKet26jq16GzhOKFo8hEh9H/AXoXwqkfieFUrBKTaXmdC8ULR5LoC6MDdvpkcZrH2wAnqOrdeQ+3aY4Tqnm4OXYB7LdEQ6qr6jP5yJWOxGgqVojXtLyZisaxOlNVbwicI6/Tm+ZadKSRPNXtDOz3CP4lKjbqcA/s72mMqs4SkR8D17NxZ+cWVZ2YkxcM3SmlEG7AltiQ5JXRbR3WA7tlynrBO1JFOeI0O9li4JSk+9Owns5bJS1riV0knVgA7+0L2NSr22Dntm8B5mO9srdI+rw8i7WyCv75zeB3ylnHrzpk2BVrDZj4XH4SfaF9ihXnf2PDqS8CdiyA9+x1YM+k+62ijBuinCvY2IFxq1A5o2zdsR2vL6P3dQXWUnFplO/W6O+uEpsCue6vGfo/qBBuWOuZZVhz2wOAq6IP8Gxgt6T1CqVwTMdanlyEtc1Ovv0o+kCflVgWOGtspmKNcsRpWt6dMrz1D/25JcD0pnXMG5upbrGRAx7CZioEazCxFBiVst49wLScvGbo/6BCuGHNb69MWbY98DLwNdA1WlYohUOwebq/woZB+EHSYy2jD/1mYxYFyvoWcG3S/c+Bs9Osdz7wdQHkXZLy5bBd9H4em7LeCQVQOBJHlzXdCuHIcwFwVtL9TlGu01PWuwj4qAA+B6mF42ugLM16wYfKwZoLH5N0v1WU/+iU9X6ONfap82t6Pw7TkZSOfar6pYj8DKvSk0SkN4XRvhy1T8FIEXkQ+BPwbjRJy5/CJkvrr8A4EfkcuJuNU7EuwU5PJToAFsJUrLBxWt5XsaOl5Gl5X1TrbFko0/J+C7wI3FnDeodjTc1Dyvv0pjlWyFPdrmbTjqCJf6cO39MM6zRaZ144zAJgN+wI43tqbbPPFpEh2KFg0IviqVR1GXCliIzE2pbPBv5GAc2JrPGaihXiNS3vW9h1uKerWyma+yS0vE9vmgNxmer2NeCPIjI7yvIPbFTq30SjYXwbjWf3a6zQ1Zm3quL7QcJ2VtWjqlnnd9jesmrAweKqIyJnY9Nb7ogNZhZ8WssEiclUrBCraXn/APRT1dQOq6nr/RS4XlW75SdZ2gx5n960LiRGU92KyK7Y8DiJz8Ec7Cj+Yawn/lygM7Yj1E1V367za3rh+L4pW0/gL1rNVLEicg52rvuiqtYJLTqN0hxYqaqVofM4B2GmN80HKZCpbqP+W4dhLQEnq+rqqCPzJWzc2blPVefn5PW8cDjnnMtGIYxC6eqJiNwhIqNC58hEnLJC/PI6l0t+cTwLInIH0EBVLw6dJUPdiM/OQZyyQozyisgk7OzCz0JnqUmcskK88uYyqxeO7MTmywJAVXcNnSFTccoKscsrxOdzG6esEK+8Ocvq1ziKWNQMs62qhp5opkZxygrxy+tcLsWlUhYEEWkSzXERFydiM5nFQZyyQozyisgWcfncxikrxCtvLrN64chObL4sXGkQkStE5BMRWS0i74jIeWlW+wkF8LmNU1aIV958Z/VrHDEkIpm2GU/X0zWv4pQV4pU36vA5DJvW9r/Y9KFjROQU4FxVzcnwErkQp6wQr7whsnrhIF5fFpGfYkN3vF/DeoUw1EScskK88v4S+Ieqfj8OVTS+2jhgioicpKpLgqXbVJyyQrzy5j2rXxwHRGQ9mX1Z7AAcHHrIERF5B5ilqj1rWO8M4IGQeeOUNcoRm7wi8i3wC1WdmrK8MzZfSENsPK3tgNc9a+bilDdEVr/GYWYC76nqmdXdgH+FDhqZBhySwXrJgwiGEqesEK+8y7FB9jahqnOw0xWLgTeAA/MbK604ZYV45c17Vj/i4PuByo5T1U41rNcDm4M6aMEVkV2AvVX1iRrWa4o1GU0dujpv4pQ1yhGbvCLyOPCtqp5bxeNNsYHujifw4JxxyhrliU3eEFm9cBCvLwvnEkTkTGAQcFJVg3OKSEPgdmxwzh/kM19KjthkjbLEJm+IrF44nHPOZcWvcTjnnMuKFw7nnHNZ8cLhSoqIXCgiM0TkWxFZKiL/FZF6aS0nIj8UketEZJsM1r1ORDTptkBEJkTX32ra9sJomxa5Se5c9bxwuJIhNv3vncBE4HTgfOBx4OR6eskfAtcCNRaOyHKga3T7JbAfMFlEmtew3dPRNhW1zOlcVrznuCslVwIjVPX3ScueFJHrQwVKsV5Vp0X/niYi84BXgBOAh1JXjlrKNFTVr4Gv8xfTlTo/4nClZBvgy9SFmtS0UEQ6R6d9zhGRe6JTWl+JyLWp24nI0SLypoh8JyKLROS2xOkiETkKeDJa9bPoOedkmXdG9LNz9JxjRGS6iJwqIjOB74CD052qEpGmInKziMwVkTUi8pmI/CUl/yUiMjN6fK6I/BrnMuBHHK6U/AcYEO3JP1XD+D1/B54CzsDGr7pWRBar6q0AIrI38BzwAtAD6Aj8FdgZG97hP0RjCGGnxRYCa7LM2zn6+WXKspuBG6LlnwGbXAcREcFOwXUFbsQK0A7AEUnr/Ar4c/RcU4EuwI0iUqGqt2SZ05UaVfWb30riBvwI+BQbLmQDNtTMDcDWSet0jh5/PmXbO4AvsKmDAcYDs7FTRYl1zoq27RrdPym63zmDbNdhQ0M0im4/BKYAK4D20TpjoufbL2XbC6PlLaL73aP7J1fxWlsDK4FrU5YnilHDmvL6rbRvfqrKlQxV/R+wJ3Yx/DZsrKk/ANPTtEh6NOX+I0AHYMfo/kHAo6pambTOBGA9cHgtI7YG1kW3D7Gjl56qujBpnS9U9e0anudo4ButeiSErkBz4CERaZS4AS8C7dj4OzqXlp+qciVFVddg1x6eBBCRi7GWVhcDQ5JW/Spl08T99sC86OeilOeuFJElwLa1jLccOAY7WvgSWKCqqUM7LNpsq821xk6NVSUxIN7MKh7vCPiwOq5KXjhcSVPVUSJyM7BHykNtq7i/MOnnJutErZxaA2nHC8rAelWdXsM6mYwRtAQrbFVJ5DuJ9IXowwxew5UwP1XlSoaIpBYDRGQ7oCWbf4GelnI/cYF7fnT/TeC0qFgkr9MIeDW6vzb6me9JnyYD24rISVU8/gawGuigqtPT3L7NX1QXR37E4UrJu9EQ1M9jp546YS2fKoCxKevuHQ23PwFrVXUxMFBVN0SP/wmbpvMxEbkduy7wN2Ciqr4RrZPYc79URMYDFar6bv38apt4AevkeJ+I3IC18GoP/FRVL1XVZSJyHTBERDoBL2M7kT8EuqlqatF0bhNeOFwpuQE4BRiKXYf4EngduwD9Wcq6v8ZO5UzA+kvcCHzfTFVVZ4rI8ViT1kew1k/3R9sl1pkrIr8ErgIGYEcrnevjF0umqioip0WZy7CZ3xYA9yWtc7OILMCG474a+x0/Ah6o73wu/nxYdeeSRNNtfoZNxflU2DTOFSa/xuGccy4rXjicc85lxU9VOeecy4ofcTjnnMuKFw7nnHNZ8cLhnHMuK144nHPOZcULh3POuaz8f8yImBZlofgsAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot exact payoff function (evaluated on the grid of the uncertainty model)\n",
    "x = uncertainty_model.values\n",
    "def payoff(x):\n",
    "    if x <= strike_price_1:\n",
    "        return 1\n",
    "    elif x < strike_price_2:\n",
    "        return 1+strike_price_1 - x\n",
    "    elif x < strike_price_3:\n",
    "        return 1+x - 2*strike_price_2 + strike_price_1\n",
    "    else:\n",
    "        return 1\n",
    "y = [payoff(x_) for x_ in x]\n",
    "plt.plot(x, y, 'ro-')\n",
    "plt.grid()\n",
    "plt.title('Payoff Function', size=15)\n",
    "plt.xlabel('Spot Price', size=15)\n",
    "plt.ylabel('Payoff', size=15)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "exact expected value:\t0.7402\n"
     ]
    }
   ],
   "source": [
    "# evaluate exact expected value (normalized to the [0, 1] interval)\n",
    "exact_value = np.dot(uncertainty_model.probabilities, y)\n",
    "print('exact expected value:\\t%.4f' % exact_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate Expected Payoff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set number of evaluation qubits (=log(samples))\n",
    "m = 6\n",
    "\n",
    "# construct amplitude estimation \n",
    "ae = AmplitudeEstimation(m, butterfly)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# result = ae.run(quantum_instance=BasicAer.get_backend('qasm_simulator'), shots=100)\n",
    "result = ae.run(quantum_instance=BasicAer.get_backend('statevector_simulator'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact value:    \t0.7402\n",
      "Estimated value:\t0.6413\n",
      "Probability:    \t0.4953\n"
     ]
    }
   ],
   "source": [
    "print('Exact value:    \\t%.4f' % exact_value)\n",
    "print('Estimated value:\\t%.4f' % result['estimation'])\n",
    "print('Probability:    \\t%.4f' % result['max_probability'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEPCAYAAAB/WNKuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAdbklEQVR4nO3de7hcVZ3m8e/LRUgaiEcu4SJNANE0GekLAck0QpB76BEMYDLo9AMNBhwUuh+kuYgQ0OERHAF9kAd4tMMwaiItNNNIQgiQEwj3IEGQBAhyB0G6D8QQiEB+88fagZ196pyz67LrpM55P89TT1WtvfaqtVcq9Ttrr73WVkRgZmbWausNdgXMzGxocoAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4xZPyRNl9Sdve6WNL3O/SdKimJZfeS9SdKj/Wy/XNIbkjYq+dmfkBSSDqmnzmat4gBjtu6YCfwXSbsWN0haHzgKuCEiVrW9ZmYNcIAxW3f8P2Al8N9rbNsPGE0KQmYdwQHGrEGSJkj6d0mvSHpL0mJJX2q0vIh4C7gJmFJj81TgNeCO7LO3kzRD0jOS3pb0pKTzJW3YT303yE6ZnVRI/46k3xfSdpD0C0k9klZKmiNpl0aPzYanDQa7AmbrsoiYnns9sbB5B+Bu4ErgHeBvgRmSVkfEzGyfbkDFsvoxE5giafeIeAggCxqTgZ9FxPtZvi2B14F/BN4AxgLnAVsAJ9d5mGuRtEV2XK8C07JjOxuYJ+lTPkVnZTnAmDUoImateS1JwJ3Ax4Gv0PiprDmkgDEVeChLOxjoypcZEYuBxbnPvxt4G7hS0qkR8V6Dnw9wGrARsH9EvJGVfw/wLHAscFUTZdsw4lNkZg2S1CXph5KeA97NHtOATzZaZkT8CbgB+GIWtCCdMnsOuDf32etJOk3SEklvZ5/9f4ARpCDXjAOAucCK7LTaBsCbwK+B8U2WbcOIA4xZ464h/fh/DzgI2AP4F2DjJsudCfw5MEHSxsDhwKxYe+nz04CLgH8FPg/sCZySbWv287cAvsSHQXPNYx9g+ybLtmHEp8jMGpD98P8dcHJEXJlLb8UfbfNJ4x9TgW2ATel9yu1oUtA5N/fZuw1Q7vvAe8BHCuldhff/CTwMXFijjOUDfIbZBxxgzBqzEekMwAcD3pI2JfUmmrrJUkS8L+k6UhDZDlgSEY8Uso3If3am3yvYIiIkvQT8Ra7O6wP7F7LeTuo1PeoBfWuGA4xZAyLiTUkPAudKWg6sBs4kjVVs1oKPmAl8HfgC6eqwonnAVyUtAn4H/D0wpkS5/wZMk/QIaVznK8DIQp7/DRwD3CHpcuBlYGtgX6A7Iq6r+2hsWPIYjFnjjiH9uF8L/AC4PnvdtIi4l3TVlqh9Rdp5wHWk01gzgbeAfypR9LmkiwguBGYAiyjUOSJeA/YClgGXAbeSxns2BfpcysasSO2+ZbKkTwCnAxOAccBdNeYX1NpvFOnLfgQpMP4KOCUi/qOQ73DgO8AupP/850fEL1p5DGZmNrDB6MGMAyYBTwBP1rHfdcBE4ATStfh7ADfmM0jam/RX5HzgUOBmYKakg5qttJmZ1WcwejDrRcTq7PUvgS0G6sFImgDcA+wbEXdmaXsC9wMHRsRtWdpcYMOI+Fxu39nAZhGxdxXHY2ZmtbW9B7MmuNTpUODVNcElK+cB4JlsG9kS5vuRejp5s0jzCUY1VmMzM2tEpwzyjwWW1khfkm0D2BnYsEa+JaTjbHh2tZmZ1a9TLlPuIq3PVNQD7JTLQ418PYXta5E0jbS8ByNGjNh9++2bm6i8evVq1luvU+L24HJblee2Ks9tVV4r2urJJ598PSK2rLWtUwJMZSLiauBqgPHjx8eiRYuaKq+7u5uJEye2oGZDn9uqPLdVeW6r8lrRVtlafDV1SpjvAWqNoXTxYQ9lzXMxX1dhu5mZtUGnBJilfDjWkpcfm3matCBfMd9Y0izrei6JNjOzJnVKgJkDbJ3NcwFA0njS+MscgGzNpPmk9ZvypgD3RsSbbaqrmZkxCGMwkkaSJlpCWshvM0lHZe9nR8RKScuABRFxPKRlMyTdClwr6RukHslFwMI1c2Ay3wa6JV1GmoQ5KXscUvmBmZnZWgZjkH8r0j0s8ta835G0/tIGwPqFPFOAS0n32/hgqZh8hohYmAWr7wBfJc2TOSYibm1h/c3MrIS2B5iIeJbsHuX95BlTI+0N4Ljs0d++N1JYQsbMzNqvU8ZgzMyswzjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSXaHmAk7SrpdkkrJb0s6QJJ6w+wz3RJ0cfjrFy+a/rIM7b6IzMzs7wN2vlhkrqA24DHgcOBnYHvkwLdOf3s+mPglkLaEcAZwJxC+lLguELas43V2MzMGtXWAAOcBIwAJkfEcmCepM2A6ZIuztJ6iYgXgRfzaZK+BSyNiMWF7G9FxH0V1N3MzOrQ7lNkhwJzC4FkFino7Fu2EEmbAwcCM1tbPTMza5V2B5ixpFNYH4iI54GV2bayjgQ2pHaA2VXSckmrJC2UVDpwmZlZ67T7FFkX8EaN9J5sW1lTgV9HxFOF9IeB+0ljPFsCp5FOw+0dEQ/UKkjSNGAawOjRo+nu7q6jGr2tWLGi6TKGC7dVeW6r8txW5VXdVu0OME2TtA3pdNoZxW0R8YNC3tnAb4GzSRcF9BIRVwNXA4wfPz4mTpzYVP26u7tptozhwm1VntuqPLdVeVW3VbtPkfUAo2qkd2XbyvgiIOAXA2WMiJXAbOBvylbQzMxao90BZimFsRZJ2wMjKYzN9GMqsDAiXiiZP7KHmZm1UbsDzBzgYEmb5tKmAG8DCwbaWdIYYC9KXj0maQRwGPBQvRU1M7PmtDvAXAmsAm6QdEA2wD4duCR/6bKkZZJ+UmP/qcB7wL8WN0gaJekuSSdK2l/SFGA+sC1wYQXHYmZm/WjrIH9E9EjaH7gcuIl0RdmlpCBTrFet5WOmArdHxOs1tq0C/kBaEWAr4B3gXmDfiFjUkgMwM7PS2n4VWUQ8DnxugDxj+kj/q372eQeY3FTlzMysZbyaspmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVoq4AI6nW8i1mZma91NuDeUnSxZL+opLamJnZkFFvgLkSOAp4TNL9kqZJ2qyCepmZWYerK8BExPSI2Ak4EHgCuAR4RdLPJB1QRQXNzKwzNTTIHxF3RMTfA1sDXwc+BcyV9Kyk6ZK2bWUlzcys8zR7Fdl4YB/SbZB7gLuAE4Blkr7cZNlmZtbB6g4wknaQdJ6kp4HbgW2AfwC2jYj/AewAXAV8r6U1NTOzjlLXDcckzQc+C7wEzABmRMRz+TwR8b6knwOntqyWZmbWceq9o+VrwCRgXkREP/kWAzs2XCszM+t49Z4i+xFwT63gImkTSfsARMS7xZ6NmZkNL/UGmPnArn1s+1S23czMrO4Ao362bQKsbKIuZmY2hAw4BpOd9pqYSzpB0iGFbBsDhwGPtq5qZmbWycoM8n+GNJkSIICjgfcKef4ELAVOb13VzMyskw0YYCLie2RzWiQ9A3whIhZXXTEzM+tsdV2mHBG+9NjMzEopMwYzCVgYEcuz1/2KiNktqZmZmXW0Mj2YXwF7AQ9kr4O+ryYLwDclMzOzUgFmR+CV3GszM7MBlRnkf67WazMzs/6UGYMZWU+BEeHJlmZmVuoU2QrS2EpZHoMxM7NSAeYfqC/AmJmZlRqDuaYN9TAzsyGm2Vsmm5mZ1VRmkP8B4NiIeFzSgwxwuiwi9mxV5czMrHOVGYP5LfB27rXHY8zMbEBlxmCOy70+ttLamJnZkNHwGIySLSX1dxMyMzMbpuoOMJImSboHeAf4PfCOpHskHdby2pmZWceqK8BIOhG4iTT58lTSzcdOzd7/e7bdzMysvvvBAGcDV0XE/yykXynpSuCbwFUtqZmZmXW0ek+RbQ78Wx/brgc+NlABknaVdLuklZJelnSBpH6Xl5E0RlLUeMyqkfdwSY9KekfS45KmlDoyMzNrqXp7MPOBfYF5NbbtC9zZ386SuoDbgMeBw4Gdge+TAt05JT7/G8DdufevF8rfmxTorgBOASYBMyX1RMStJco3M7MWKTPRctfc2x8CP5a0OXAj8BqwFfAF4FDghAGKOwkYAUyOiOXAPEmbAdMlXZyl9eeJiLivn+3fAu6MiFOy9/MljQPOBRxgzMzaqEwP5jHWnlwp4MTsUby75S30v5ryocDcQiCZBVxE6gHdVKI+NUnaCNiP1HPJmwXMkDQqIt5stHwzM6tPmQCzXws/byxwRz4hIp6XtDLbNlCAmSHpY6Se00zgmxGxZpWBnYENgaWFfZaQTsF9EniwueqbmVlZZWbyL2jh53UBb9RI78m29WUV8CPSaa7lwETgDFJQOTxXNjXK7ylsX4ukacA0gNGjR9Pd3d1f/Qe0YsWKpssYLtxW5bmtynNblVd1W9U7yP8BSesBGxfTq7ijZUS8Anwtl9Qt6VXgCkl/GRGPNFH21cDVAOPHj4+JEyc2Vdfu7m6aLWO4cFuV57YqbzDaasyZN6/1/tnvdsa886rbqt6JlpJ0hqRlwLvAH2s8+tMDjKqR3sWHPY2yfpk9754rmxrldxW2m5lZG9Q7D+YU4EzgJ6TB/f8FXAA8CTxLdqqpH0tJYy0fkLQ9MJLeYycDicLz06SgN7aQbyywOqujmZm1Sb0B5ivAecDF2fsbI+J8YBwpQOwywP5zgIMlbZpLm0K6HUC9Yz1HZc8PAUTEKtI8naML+aYA9/oKMjOz9qp3DGZHYHFEvC/pXeCjABGxWtIVwI9JPZy+XEnqBd0g6SJgJ2A6cEn+0uXsFNyCiDg+ez8d2JQ0yXI5sA9wOnBDRPwmV/63SeMzl5Hm6UzKHofUeZxmZtakensw/wFskr1+Hvjr3LYu0iTKPkVED7A/aa7MTcD5wKWkXlHeBqw9n2YpaZ7MDGA2cAzwvew5X/5CUs/mAGAu8HngGM/iNzNrv3p7MHcDe5B+5H9OmoH/MeBPwMnA7QMVEBGPA58bIM+YwvtZpAmTA4qIG0m9FzMzG0T1BpjpwHbZ6wtJp8iOJfVc5gFfb1XFzMyss9UVYCLiCeCJ7PUq0r1gTq2gXmZm1uGamWj5cWAb4OWIeKl1VTIzs6GgkVsmf1XSC8BzwP3A85JelFS8CZmZmQ1j9c7kPxe4nDSf5TBgfPY8B/hhtt3MzKzuU2QnAxdGxLcK6bdka4OdTJrZb2Zmw1y9p8hG0PddKxdQY/FLMzMbnuoNMDcCk/vYdiTwq+aqY2ZmQ0WZWyZPyr2dA1wsaQy9b5k8Dvjn1lfRzMw6UZkxmF/R+9bI2wEH18j7U9KdJs3MbJgrE2B2rLwWZmY25JS5ZfJz7aiImZkNLXXP5Je0AWlAf2/gY8B/AneRls5/r7XVMzOzTlVXgJG0FXArsBvpDpavAhNI818ekXRQRPyh1ZU0M7POU+9lypcAmwN7RcROETEhInYCPpOlX9LqCpqZWWeqN8BMAs6IiAfyiRHxIHAWadkYMzOzugPMRsAf+9j2R+AjzVXHzMyGinoDzH3AGZL+LJ+YvT8j225mZlb3VWSnAfOBFyTdShrk34o06VLAxJbWzszMOlZdPZiIWAzsAlwNbAkcSAowVwK7RMQjLa+hmZl1pNI9GEkbAnsCz0TEmdVVyczMhoJ6ejDvA3cAYyuqi5mZDSGlA0xErAaeAraurjpmZjZU1HsV2TeBcyV9uorKmJnZ0FHvVWTnkGbsL5b0EukqsshniIg9W1Q3MzPrYPUGmMeyh5mZWb9KBRhJI0jLxDwG/B64LSJerbJiZmbW2crcMnkn4DZgTC55uaQvRsStVVXMzMw6W5lB/ouB1cBngZHAOOBh4KoK62VmZh2uTICZAJwTEXdHxDsRsQQ4EfhzSdtUWz0zM+tUZQLMNsDvCmlPk9Ye85wYMzOrqexVZDFwFjMz68+YM29e6/2z3x3at9AqG2DmSnqvRvrtxfSI2Kr5apmZWacrE2DOr7wWZmY25AwYYCLCAcbMzOpW71pkZmZmpTjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVou0BRtKukm6XtFLSy5IukLT+APvsIWmGpGXZfk9IOk/SxoV80yVFjcch1R6VmZkV1XvDsaZI6iIt/f84cDiwM/B9UqA7p59dp2R5LwKeAnYDvp09H1nI+yZQDChLmq27mZnVp60BBjgJGAFMjojlwDxJmwHTJV2cpdXy3Yh4Pfe+W9I7wFWSdoiI53Lb3ouI+6qpvpmZldXuU2SHAnMLgWQWKejs29dOheCyxsPZ87atq56ZmbVKuwPMWGBpPiEingdWZtvqMYF0I7SnC+kflfS6pHclPSxpcsO1NTOzhimifSvxS3oXOD0iLiukvwhcGxFnlyxna+A3wOyIODaX/mVgK1LvZlPSjdEmAUdGxA19lDUNmAYwevTo3WfNmlXvYa1lxYoVbLLJJk2VMVy4rcpzW5U3GG316EtvrvX+09uNaipfu7Sirfbbb7+HImJ8rW0dF2AkfYR0ocDHgd0joqefvALuAUZExF8NVPb48eNj0aJFA2XrV3d3NxMnTmyqjOHCbVWe26q8wWirsvd5WdfuB9OKtpLUZ4Bp9ymyHqBWyO7KtvUrCxjXAuOASf0FF4BI0fMGYLeBLoU2M7PWavdVZEspjLVI2h4YSWFspg+XkS5vPjAiyuSHdDdO35HTzKzN2t2DmQMcLGnTXNoU4G1gQX87SjoL+Brw5YhYWObDsh7PkcAjEfF+Y1U2M7NGtLsHcyVwCnCDpIuAnYDpwCX5S5clLQMWRMTx2ftjgAuBa4CXJO2VK/PpiPhDlm8BcD2pN/RnwFeAzwBHVHtYZmZW1NYAExE9kvYHLgduAt4ALiUFmWK98mMmB2XPx2aPvONIgQdgGfCPwDakS5h/DRwWEXNaUX8zMyuv3T0YIuJx4HMD5BlTeH8svQNLrf2Ob6JqZmbWQl5N2czMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklNhjsCpiZrYvGnHnzWu+f/e5hg1STzuUejJmZVcI9GDOzdUSx1wSd3XNyD8bMzCrhAGNmZpVwgDEzs0q0PcBI2lXS7ZJWSnpZ0gWS1i+x3yhJMyT1SHpT0s8kbV4j3+GSHpX0jqTHJU2p5kjMzKw/bR3kl9QF3AY8DhwO7Ax8nxTozhlg9+uATwInAKuBi4Abgc/myt8buB64AjgFmATMlNQTEbe29GDMrKMMhQH0Trt0ut1XkZ0EjAAmR8RyYJ6kzYDpki7O0nqRNAE4CNg3Iu7M0l4C7pd0QETclmX9FnBnRJySvZ8vaRxwLuAAY2bWRu0OMIcCcwuBZBapN7IvcFM/+726JrgARMQDkp7Jtt0maSNgP1LPJW8WMEPSqIh4s0XHYWaDaCj0RlppXW2PdgeYscAd+YSIeF7SymxbXwFmLLC0RvqSbBuk020b1si3hHQK7pPAg41V26xxzfznX7PvaZ9+j2PPvLnu/Zr5zGb3K7vvuvrj2OlqtWtRvd+reikiKim45odJ7wKnR8RlhfQXgWsj4uw+9psHvBURRxTSfwrsFBH/VdLfAguBv46Ixbk8nwCeAg6uNQ4jaRowLXv7KeCJhg8w2QJ4vckyhgu3VXluq/LcVuW1oq12iIgta20Y9jP5I+Jq4OpWlSdpUUSMb1V5Q5nbqjy3VXluq/Kqbqt2X6bcA4yqkd6VbWtmvzXPxXxdhe1mZtYG7Q4wS/lwzAQASdsDI6k9xtLnfpn82MzTwLs18o0lXdb8ZAP1NTOzBrU7wMwBDpa0aS5tCvA2sGCA/bbO5rkAIGk8sFO2jYhYBcwHji7sOwW4t41XkLXsdNsw4LYqz21VntuqvErbqt2D/F2kSZaPkS5N3gm4BLgsIs7J5VsGLIiI43Npc4FdgG/w4UTL1yKiONGyG7icNAlzUpb/EE+0NDNrr7b2YCKiB9gfWJ90SfL5wKXAeYWsG2R58qaQejn/AlwLPAR8oVD+QuAo4ABgLvB54BgHFzOz9mtrD8bMzIYPr6ZcgxfkLK+RtpK0R9ZOy7L9npB0nqSNC/mmS4oaj0OqPapqNNhWY/pog1k18g7371Vf35eQdFYu3zV95Kl1IdE6T9InJF0l6TeS3pfUXXK/yn+vhv08mCIvyFleE201Jct7EWkS7G7At7PnIwt53wSKAWVJs3Vvtya/V5DGEu/OvV9rcpy/VwD8GLilkHYEcAbZxUA5S4HjCmnPNlbjQTeO9O99H2k1k7Kq/72KCD9yD+As0pyZzXJp/wyszKfV2G8CEMA+ubQ9s7QDcmlzgTsK+84GFg72sbexrbaokTYta6sdcmnTgdcH+zgHua3GZO3ydwOUP+y/V32UdTOwpJB2DbBosI+zhe21Xu71L4HuEvu05ffKp8h662tBzhGkBTn726/XgpzAmgU5yS3IeV1h31nABEm1JpOuyxpqq4iotTTFw9nztq2r3jql0e/VgPy9qi073XMgMLO11Vu3RMTqBnZry++VA0xvvRbWjIjnSX899XeOtlULcnaSRtuqlgmkbvrThfSPSnpd0ruSHpY0ueHaDq5m22pGdn79FUmXSBqR2+bvVW1HktqlVoDZVdJySaskLZTUVJDvQG35vXKA6a0LeKNGeg8fLjvT6H5rnov5egrbO0WjbbUWSVuTzq3/34h4LbdpGenUyNGkH4uXges7NMg02largB8Bx5Mu8b8K+Crpr8h82dQof1h/r4CpwK8j4qlC+sPAacB/A75EmhIxT9KeDdS1U7Xl98qD/DaoJH2E1AVfAfxTfltE/LSQ9ybgHtIN5G5oVx0HU0S8Anwtl9Qt6VXgCkl/GRGPDFLV1mmStiGdTjujuC0iflDIOxv4LXA26aIAaxH3YHrzgpzlNdpWAEgSadLsOGBSpIm4fYo0wngDsFuZy8bXMU21VcEvs+fdc2VTo/xh+b3KfBEQ8IuBMkbEStLA9d+UreAQ0JbfKweY3rwgZ3mNttUal5EuQz08Isrkh3SVSyfODm62rfKi8OzvVW9TSVc6vVAyf6d+rxrVlt8rB5jehsOCnK3SaFuRTXz7GvDlSEv8DCjr8RwJPBIR7zdW5UHTcFvVcFT2/BD4e1UkaQywFyWvHssumDiMrD2Hifb8Xg32Ndzr2oPU/XsFmEda02waaXzgO4V8y4CfFNLmAr8DJpPO5T4B3FXIszfwHumv94nAxaS/Bg4a7GNvV1sBx5D+WpxB+iHIP7bM5VtAmtx1EGndudlZW31+sI+9jW01nTTJcHK23wWkH9rr/b3q/X8wSz+T9Jd3rflWo4C7gBNJF01MIU1QXAWMH+xjb7C9RpL+6DgKuJc0nrTm/ci+2qodv1eD3jjr4gPYFbgj+4/8CmmW+fqFPM8C1xTSPpr9aL4BLAd+3seX/AjSitKrSN3RqYN9zO1sK9JEt+jjcWwu30+y/wBvA29lPwyHDvYxt7mtpgKLSCsa/Cn7obgA2Mjfq97/B7P0xcAtfZS7MWkc74Wsnd4kzf7fa7CPuYm2GtPP/6cxfbVVO36vvNilmZlVwmMwZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwS/x+pgYWP5m0jEgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEPCAYAAAB/WNKuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3debhcVZnv8e+PoSUQCGFIAEUitpgL6lUSGbppSQRUgm0QxNBoeyOQxG5stB+hGUQJolwBGbq1vSSoIK1NHFC6lSEESII4QYLQKCQYJIwyB0JIgATe+8faBzY7VXWq6tSuOvvw+zzPfurU2kO9q8459dZea6+1FRGYmZl12ga9DsDMzIYmJxgzMyuFE4yZmZXCCcbMzErhBGNmZqVwgjEzs1I4wVhHSZopKeosH2/yGLtkx9myUD41O87wcqJvLo4BHvPHkhY0sd1Gkj4r6TZJayStkHSlpH3afN3B8p5OLfxNPCpprqTdm9h3QrbP27oRqw2cE4yV4Wlg7xrL1U3uvwtwKlD8YL8iO87qzoTZdhylkrQhcDlwBvDfwCRgKvAisEDSEW0cdrC8p33em73uDGBbYL6kHfrZ55Zsn7tLjs06ZKNeB2BD0rqI+E2nDxoRjwGPdfq4g9A/AQcBB0ZEPin/l6Q5wGxJCyPiwYG+UA/f05sjYhWApEXAvcDHgLOLG0oS8LqIWAl0/O/KyuMzGOsJSSdJWibpOUmPSLpa0naSJgA/yza7J2sSWZ7t86rmHEljsueHS7pI0kpJD/Q1xUn6F0kPSXpM0pmSNsi9/lhJcyTdL2m1pD9kTVIbZOvrxpGtf2O2/5PZ/nMlvbVQxx2zZq01kpZLOrrJt+czwPxCcunzeWAT4Kjc6yyX9DVJX5D0sKRVkr4vaUR/danVRCZpG0nflfREVrcFksYX6tb3mv+cvecrsvej5bO9iLiflOTGZMeeKelxSftIuhl4DjisVhOZpA2zv6W7JD2fxXJxIdbJkhZlf2sPSzpL0satxmmt8xmMlULSen9bEbEuW/cJ4GTgBOAPwNakJpPNSM0gxwFfAw4B/gw838/LnQl8HzgUOBL4rqR3ATtlz8cBXwZ+B8zJ9nk9sDTb7xngncBpwDDg/zaKQ9JWwI3AE8CnSM1LJwLXStolItZk37r/C9iGlAyey46/FfDHBu/bjqQP2vNqrY+IuyXdDrynsOrvgGXANGB74CzgW8BhjepSx+XAX2b7PA4cT2rCeldELMtt91Hgf4DpwBuAc0nNev/Y4NjrkbQ56X15OFe8KfDdrB53AQ9l9SqaBXwi225hdpxDc8f+KHBptt3JwJtJv98NsvpZmSLCi5eOLcBMIOosY7JtvgFc1uAYH8xvnyufmpUPz56PyZ5flNtmC2At6UN8w1z5TcAP6ryeSF+2Tgb+1EQcp5OSy1a5spGkvqdjsueTsn33zG2zE7AOWNCg7ntl+01usM3lwJ2558uBJ/vel6zsY8BLwP9q8T39QPZ839w2m5HOMGYVXvNuYKNc2fnAw/38ffS93ojsPd8R+EH2vryz8Dc0ubDvhKz8bdnzsdnzYxv8Xu/N/31k5UcCa4Cte/3/MtQXn8FYGZ4G9q9R/lD2eCtwlKTTSJ3MiyPixQG83nV9P0TESkmPAQsLx1wGvLHviaRNgJNIH8RvBDbOrdsosrOtOvYH5gErc2dqzwCLgb6mpD2ARyLit7nY7pW0uI36NWNeZH0amZ+SPmDfDdzZwnH2AB6NiIV9BRHxrKSfA8Ur2OYX3qc7gFGSNo6Itf28zlO5nx8HjoyIW3NlAVzVzzEmZo8X11m/C+l3+8PCGfX1pGbGt5HOeqwkTjBWhnURsajB+u8Am5OaVr4IPCHpAuDUNhPNU4XnL9Qp2yT3/EzgaFKz1S3Z9pOBU7LtVlHfNqQzjSk11vUlu+2AR2usf5RU93r6Ou53arDNTrnt8sd9WUSslrSK2s1KjWxfPFbmEVLzU16t91jA60hnkY28h9S0+Dhwf0S8VFi/IiJe6OcYWwPPRur8r2Wb7PHKOut37Of4NkBOMNZ12YfJecB5WZ/Dx4CvAA8AF3QpjMOAr0fEWX0Fkg5qct8nSZcPn15j3TPZ48PAqBrrR5GaZ2qKiPuzDvgPAf9WXC/pTaRv3sXXHlXYblNgOKm/pRV/Lh4rM5pU7075XeGMq6iZ+4g8AWwmaYs6SaYv3umk/reie5p4DRsAX0VmPRUR90fEV0lNWLtmxX3fXDepvVdHDCPX0a009uTwwjb14rgO2A34Q0QsKixLs21uBkZL2jP3Gm8E+h1QCPwrsJ+k99VY95Us7m8Xyg/QqwdLfpj0Id13Jtnse/pbUjPXyxcRZMnqINKFDYPJ9dnjJ+qsX0o60xtT4/e0KCKe6E6Yr10+g7EybCRprxrl90fEg5Jmkb5d/obUXzMReAvpqjJIHwwAM5TGfayOiNs7HOM84BhJy7JYjiE17eTVi+Nc4OPA9ZK+TvoQGw3sC9wYEZeSmmVuA34k6QRSUjiN2s1PRV8n9fP8VNLXgAWkZrWjSJ31fx/rj4FZA1wh6WxSM9fZwE8j4o5+6vIqETFX0q+AH0g6kXSWcBwpIa83RqWXImKppNnAOZJGATeQBpJ+JCIOj4iXJH0O+A9JW5D6dF4AdgYOzrbr9gDT15ZeX2XgZWgtNL6K7JRsm6nAL0kf7KtJl7oeVTjO50hXAK0Dluf2q3UV2QcL+y4HvlYouxhYlHs+mtQRvpLUv3AW6RLfl49fL46sfAfgomzf57PX/B6wW26bN5JmL1iTHWMG8GMaXEWW23cj4J+z92YNsIL0AblPjW2XA+dk7/0jwLOkS3O3bPU9zcq2BS7JXnMNqSP83U28x+sdq0aszWwzE3i8RvkEcleRZWUbkl39R0oeDwDfKex3IPCL7H1ZSbrI5MvkroDzUs6i7BfQNZL+knRd/d6kZoZfRMSEJvYbQboM8mBS097PSZcnPlHYbjLpj+ctpD+60yLiB52sg9lgkvXZ/DgiPK7DBpVe9MHsRhojsJQ0gKpZPyR9gzma9C3o3aTxAC9TmgjwMmA+6VvLFcClddqyzcysRL04g9kgsksSJf0Y2Ka/MxhJewO/Ig3+uiEr24PUIXlARFyblc0FNo6I9+b2vRLYIiLamoXWbLDzGYwNVl0/g4n1r3dvxoGkQWs35I5zE+kywwMBJL2O1Fn8w8K+c4C9++ZlMhtqImKMk4sNRlW5THkssKRG+Z3ZOkhzDG1cY7s7SfXcpbTozMxsPVW5THkk648ahnSVy865baix3YrC+leRNJ00EIthw4aN23HH2oN7X3rpJTbYoCr5uDHXZXAaKnUZKvWAwVWXze9KXdbP7NLed+Wy6nLXXXc9HhHb1lpXlQRTmoiYDcwGGD9+fCxaVHuGkwULFjBhwoQuRlYe12VwGip1GSr1gEFWFyk9Ll3aeLs6yqqLpHvrrRscqbl/K0izrxaN5JUzlL7H4nYjC+vNzKwLqpJglvBKX0tevm/mbtIEe8XtxpKmLW/lkmgzMxugqiSYq4DtsnEuAGR32Ns5W0dEPE8a/3JYYd8pwK8j4ukuxWpmZvSgDyabOG9S9vT1wBaSPpI9vzLSNOPLSPfzOAogIn4t6RrgEknHkc5IziTN+3Rt7vCnAwsknU8ahDkpWz5QesXMzOxVetHJPwr4UaGs7/mbSHMcbUSaYyhvCmmK9++Qmyomv0FE3Jglqy8D/0AaJ3NERFzTwfjNzLqvy4PiO6HrCSYilpNuStRomzE1yp4CPpktjfa9nMIUMmZm1n1V6YMxM7OKcYIxM6uCcePSUiGv+YGWZmaVcMstvY6gZT6DMTOzUjjBmJlZKZxgzMysFE4wZmZWCicYMzMrha8iMzOrgmnTeh1By5xgzMyqYPbsXkfQMjeRmZlZKZxgzMyqYPHitFSIm8jMzKpg/Pj0WKFZlX0GY2ZmpXCCMTOzUjjBmJlZKZxgzMysFE4wZmZWCicYMzMrhS9TNjOrgkWLeh1By5xgzMyqoGK3SwY3kZmZWUmcYMzMqmD69LRUiBOMmVkVXHhhWirECcbMzErhBGNmZqVwgjEzs1I4wZiZWSmcYMzMrBQeaGlmVgW7797rCFrmBGNmVgUVu10yuInMzMxK4gRjZmalcIIxM6sCKS0V4gRjZmalcIIxM7NSOMGYmVkpnGDMzKwUTjBmZlYKJxgzMyuFR/KbmVXBrFm9jqBlTjBmZlVQsdslQw+ayCTtKuk6SaslPSTpS5I27GefmZKiznJSbruL62wztvyamZlZXlfPYCSNBK4F7gAmA28GziElulMa7Pot4OpC2cHACcBVhfIlwCcLZcvbi9jMbJCYPTs9VuhMpttNZJ8ChgGHRMRKYJ6kLYCZks7KytYTEQ8AD+TLJH0BWBIRtxY2fzYiflNC7GZmvTNjRnqsUILpdhPZgcDcQiKZQ0o6+zZ7EElbAwcAl3Y2PDMz65RuJ5ixpCasl0XEfcDqbF2zDgU2pnaC2VXSSknPS7pRUtOJy8zMOkcR0b0Xk9YCx0fE+YXyB4BLIuLkJo9zPTAiIsYVyj8DvEDq49kW+BwwDtgnIm6qc6zpwHSA0aNHj5szZ07N11y1ahXDhw9vJrxBz3UZnIZKXYZKPWBw1WXCxIkALJg/v639y6rLxIkTF0fE+JorI6JrC7AW+GyN8geAM5o8xvbAi8BxTWy7KXAPcHkzxx43blzUM3/+/LrrqsZ1GZyGSl2GSj0iBlldIC1tKqsuwKKo85na7SayFcCIGuUjs3XN+Cgg4Af9bRgRq4ErgerdzNrMrOK6nWCWUOhrkbQj6UxjSc091nc4cGNE3N/k9pEtZmbWRd1OMFcB75e0ea5sCrAGWNjfzpLGAHvR5NVjkoYBBwGLWw3UzGxQ6Wskq5BuJ5gLgOeBn0jaP+tgnwmcG7lLlyUtk/TtGvsfDqwDflRcIWmEpF9ImiFpP0lTgPnADsAZJdTFzMwa6OpAy4hYIWk/4BvAz4CngPNISaYYV63pYw4HrouIx2usex54jDQjwCjgOeDXwL4RsagjFTAzs6Z1fbLLiLgDeG8/24ypU/7OBvs8BxwyoODMzAarcdmojMXVafH3bMpmZlVwyy29jqBlvuGYmZmVwgnGzMxK4QRjZmalcIIxM7NSOMGYmVkpfBWZmVkVTJvW6wha5gRjZlYFfbdMrhA3kZmZWSlaSjCSak3fYmZmZVu8uFKj+KH1JrIHJV0CXBQRd5YRkJmZ1TA+u2lkhWZUbrWJ7ALgI8DvJf1W0nRJW5QQl5mZVVxLCSYiZkbEzsABwFLgXODPkr4vaf8yAjQzs2pqq5M/Iq6PiE8A2wH/BLwVmCtpuaSZknboZJBmZlY9A72KbDzwHtJtkFcAvwCOBpZJ+vgAj21mZhXWcoKRtJOkUyXdDVwHbA8cCewQEX8P7ATMAs7uaKRmZlYpLV1FJmk+8DfAg8BFpKvJ7s1vExEvSvpP4DMdi9LMzCqn1cuUHwUmAfMiGl4rdyvwprajMjOzV1tUvTu/t5pg/h24pVZykTQc2D0iboiItcC96+1tZmbt6btlcoW02gczH9i1zrq3ZuvNzMxaTjBqsG44sHoAsZiZWT3Tp6elQvptIpP0HmBCruhoSR8obLYJcBBwe+dCMzOzl114YXqs0KzKzfTB7EkaTAkQwGHAusI2LwBLgOM7F5qZmVVZvwkmIs4mG9Mi6R7gwxFxa9mBmZlZtbV0FVlE+NJjMzNrSjN9MJOAGyNiZfZzQxFxZUciMzOzSmvmDObnwF7ATdnPQf2ryQLwTcnMzKypBPMm4M+5n83MrNt2373XEbSsmU7+e2v9bGZmXVSx2yVDc30wm7ZywIjwYEszM2uqiWwVqW+lWe6DMTOzphLMkbSWYMzMrNOUXVvVcCL7waWZPpiLuxCHmZkNMQO9ZbKZmVlNzXTy3wRMjYg7JN1MP81lEbFHp4IzM7PqaqYP5g/AmtzP1WkANDOznmmmD+aTuZ+nlhqNmZkNGW33wSjZVlKjm5CZmdlrVEuzKcPLk1+eAozL9l8naTHwlYi4osPxmZkZwKxZvY6gZS0lGEkzgG8C1wGfAR4FRgGHAP8t6R8jonrvgpnZYFex2yVD62cwJwOzIuIfC+UXSLoA+DzgBGNmZi33wWwN/LTOusuArfo7gKRdJV0nabWkhyR9SVLD6WUkjZEUNZY5NbadLOl2Sc9JukPSlKZqZmY2mM2enZYKafUMZj6wLzCvxrp9gRsa7SxpJHAtcAcwGXgzcA4p0Z3SxOsfB/wy9/zxwvH3ISW6bwLHApOASyWtiIhrmji+mdngNGNGeqxQU1kzAy13zT39N+BbkrYGLueVPpgPAwcCR/dzuE8Bw4BDImIlME/SFsBMSWdlZY0sjYjfNFj/BeCGiDg2ez5f0m7AFwEnGDOzLmrmDOb3vHpwpYAZ2VK8u+XVNJ5N+UBgbiGRzAHOJJ0B/ayJeGqS9DpgIunMJW8OcJGkERHxdLvHNzOz1jSTYCZ28PXGAtfnCyLiPkmrs3X9JZiLJG1FOnO6FPh8RPTNMvBmYGNgSWGfO0lNcLsANw8sfDMza1YzI/kXdvD1RgJP1Shfka2r53ng30nNXCuBCcAJpKQyOXdsahx/RWH9q0iaDkwHGD16NAsWLKgZwKpVq+quqxrXZXAaKnUZKvWAwVWXCdlju/H0oi4tD7TsI2kDYJNieRl3tIyIPwOfzhUtkPQI8E1J/zsibhvAsWcDswHGjx8fEyZMqLndggULqLeualyXwWmo1GWo1AMa12XMiY3HlS//6kElRETb720vfi8tXaacTQ9zgqRlwFrgmRpLIyuAETXKR/LKmUazfpw9jssdmxrHH1lYb2ZmXdDqOJhjgROBb5M6978CfAm4C1hO1tTUwBJSX8vLJO0IbMr6fSf9icLj3aSkN7aw3VjgpSxGM7NqiqjU3Syh9QQzDTgVOCt7fnlEnAbsRkoQb+ln/6uA90vaPFc2hXQ7gFb7ej6SPS4GiIjnSeN0DitsNwX4ta8gMzPrrlb7YN4E3BoRL0paC2wJEBEvSfom8C3SGU49F5DOgn4i6UxgZ2AmcG7+0uWsCW5hRByVPZ8JbE4aZLkSeA9wPPCTiPif3PFPJ/XPnE8apzMpWz7QYj3NzGyAWj2DeQIYnv18H/Cu3LqRpEGUdUXECmA/0liZnwGnAeeRzoryNuLV42mWkMbJXARcCRwBnJ095o9/I+nMZn9gLvAh4AiP4jezyhs3Li0V0uoZzC+Bd5M+5P+TNAJ/K+AF4BjSLMsNRcQdwHv72WZM4fkc0oDJfkXE5aSzFzOzoeOWW3odQctaTTAzgddnP59BaiKbSjpzmQf8U6cCMzOzamspwUTEUmBp9vPzpHvCfKaEuMzMrOIGMtDyDcD2wEMR8WDnQjIzs6Gg1U5+JP2DpPuBe4HfAvdJekBS8SZkZmb2GtbqSP4vAt8gjWc5CBifPV4F/Fu23szMrOUmsmOAMyLiC4Xyq7O5wY4hjew3M7NOmjat1xG0rNUEM4z6d61ciK8iMzMrR8Vulwyt98FcDhxSZ92hwM8HFo6ZmQ0VzdwyeVLu6VXAWZLGsP4tk3cD/qXzIZqZGYsXp8cKjeZvpons56x/a+TXA++vse33SHeaNDOzTho/Pj1WaEblZhLMm0qPwszMhpxmbpl8bzcCMTOzoaXlkfySNiJ16O8DbAU8CfyCNHX+us6GZ2ZmVdVSgpE0CrgGeAfpDpaPAHuTxr/cJul9EfFYp4M0M7PqafUy5XOBrYG9ImLniNg7InYG9szKz+10gGZmVk2tJphJwAkRcVO+MCJuBk4iTRtjZmbWch/M64Bn6qx7BviLgYVjZmY1LVrU6wha1mqC+Q1wgqTrI+LZvkJJmwEnZOvNzKzTKjTAsk+rCeZzwHzgfknXkDr5R5EGXQqY0NHozMysslrqg4mIW4G3ALOBbYEDSAnmAuAtEXFbxyM0MzOYPj0tFdL0GYykjYE9gHsi4sTyQjIzs/VceGF6rNCsyq2cwbwIXA+MLSkWMzMbQppOMBHxEvBHYLvywjEzs6Gi1XEwnwe+KOntZQRjZmZDR6tXkZ1CGrF/q6QHSVeRvWru6IjYo0OxmZlZhbWaYH6fLWZmZg01lWAkDSNNE/N74GHg2oh4pMzAzMwsZ/fdex1By5q5ZfLOwLXAmFzxSkkfjYhrygrMzMxy+m6ZXCHNdPKfBbwE/A2wKbAb8DtgVolxmZlZxTWTYPYGTomIX0bEcxFxJzADeKOk7csNz8zMqqqZBLM98KdC2d2kucc8JsbMrBuktFRIs1eRRf+bmJlZ0ZgTr6i7bvlXh/YttJpNMHMlratRfl2xPCJGDTwsMzOrumYSzGmlR2FmZkNOvwkmIpxgzMysZa3ORWZmZtYUJxgzMytFq3ORmZlZL8yq3th2Jxgzsyqo2O2SwU1kZmZWEicYM7MqmD07LRXS9QQjaVdJ10laLekhSV+StGE/+7xb0kWSlmX7LZV0qqRNCtvNlBQ1lg+UWyszs5LNmJGWCulqH4ykkaSp/+8AJgNvBs4hJbpTGuw6Jdv2TOCPwDuA07PHQwvbPg0UE8qdA43dzMxa0+1O/k8Bw4BDImIlME/SFsBMSWdlZbV8NSIezz1fIOk5YJaknSLi3ty6dRHxm3LCNzOzZnW7iexAYG4hkcwhJZ196+1USC59fpc97tC58MzMrFO6nWDGAkvyBRFxH7A6W9eKvUk3Qru7UL6lpMclrZX0O0mHtB2tmZm1TRHdm4lf0lrg+Ig4v1D+AHBJRJzc5HG2A/4HuDIipubKPw6MIp3dbE66Mdok4NCI+EmdY00HpgOMHj163Jw5c2q+5qpVqxg+fHgz4Q16rsvgNFTqMlTqAY3rcvuDTzfc9+2vH9Hvdn3bNGPCxIkALJg/v+l98sr6vUycOHFxRIyvta5yCUbSX5AuFHgDMC4iVjTYVsCvgGER8c7+jj1+/PhYtGhRzXULFixgwoQJ/R2iElyXwWmo1GWo1AMa16XRfV7glXu9dOx+MH03G2vzM7us34ukugmm201kK4BaKXtktq6hLGFcAuwGTGqUXAAiZc+fAO/o71JoM7NBLaLt5NIr3b6KbAmFvhZJOwKbUuibqeN80uXNB0REM9tDuhtntX4rZmZDQLfPYK4C3i9p81zZFGANsLDRjpJOAj4NfDwibmzmxbIznkOB2yLixfZCNjOzdnT7DOYC4FjgJ5LOBHYGZgLn5i9dlrQMWBgRR2XPjwDOAC4GHpS0V+6Yd0fEY9l2C4HLSGdDmwHTgD2Bg8utlplZycaNS4+LF/c2jhZ0NcFExApJ+wHfAH4GPAWcR0oyxbjyfSbvyx6nZkveJ0mJB2AZ8Flge9IlzLcAB0XEVZ2I38ysZ265pdcRtKzr0/VHxB3Ae/vZZkzh+VTWTyy19jtqAKGZmVkHeTZlMzMrhROMmZmVwgnGzMxK4QRjZmal6Honv5mZtWHatF5H0DInGDOzKqjY7ZLBTWRmZlYSJxgzsypYvLhSo/jBTWRmZtUwPpsRv0IzKvsMxszMSuEEY2ZmpXCCMTOzUjjBmJlZKZxgzMysFE4wZmZWCl+mbGZWBYsW9TqCljnBmJlVQd8tkyvETWRmZlYKJxgzsyqYPj0tFeIEY2ZWBRdemJYKcYIxM7NSOMGYmVkpnGDMzKwUTjBmZlYKJxgzMyuFB1qamVXB7rv3OoKWOcGYmVVBxW6XDG4iMzOzkjjBmJlZKZxgzMyqQEpLhTjBmJlZKZxgzMysFE4wZmZWCicYMzMrhROMmZmVwgnGzMxK4ZH8ZmaZMSde8arnn3v7OqbmypZ/9aBuh/SKWbN699ptcoIxM6uCit0uGZxgzMx6qnjWVNTTs6YBch+MmVkVzJ6dlgrxGYyZWRXMmJEeK9RU5jMYMzMrRdcTjKRdJV0nabWkhyR9SdKGTew3QtJFklZIelrS9yVtXWO7yZJul/ScpDskTSmnJmZm1khXm8gkjQSuBe4AJgNvBs4hJbpT+tn9h8AuwNHAS8CZwOXA3+SOvw9wGfBN4FhgEnCppBURcU1HK2NmldCoE70KHeh98S8vPO8zmOvQ7T6YTwHDgEMiYiUwT9IWwExJZ2Vl65G0N/A+YN+IuCErexD4raT9I+LabNMvADdExLHZ8/mSdgO+CDjBmJl1UbcTzIHA3EIimUM6G9kX+FmD/R7pSy4AEXGTpHuydddKeh0wkXTmkjcHuEjSiIh4ukP1MLMeGMqX9Lar3nsyGN6LbieYscD1+YKIuE/S6mxdvQQzFlhSo/zObB2k5raNa2x3J6kJbhfg5vbCNmtPu80zA/kg7fZrtjv63cmiXP39XvLKeq8VEaUcuOaLSWuB4yPi/EL5A8AlEXFynf3mAc9GxMGF8u8BO0fEX0n6a+BG4F0RcWtum78E/gi8v1Y/jKTpQN91f28FltYJfxvg8SaqWQWuy+A0VOoyVOoBrkszdoqIbWuteM2Pg4mI2UC/o5ckLYqI8V0IqXSuy+A0VOoyVOoBrstAdfsy5RXAiBrlI7N1A9mv77G43cjCejMz64JuJ5glvNJnAoCkHYFNqd3HUne/TL5v5m5gbY3txpIua76rjXjNzKxN3U4wVwHvl7R5rmwKsAZY2M9+22XjXACQNB7YOVtHRDwPzAcOK+w7Bfh1B64gq9YkQI25LoPTUKnLUKkHuC4D0u1O/pGkQZa/J12avDNwLnB+RJyS224ZsDAijsqVzQXeAhzHKwMtH42I4kDLBcA3SIMwJ2Xbf8ADLc3MuqurZzARsQLYD9iQdEnyacB5wKmFTTfKtsmbQjrL+Q5wCbAY+HDh+DcCHwH2B+YCHwKOcHIxM+u+rp7BmJnZa4dnU65D0jRJf8wmzVwsab8W93+XpBcl9fwa+nbqImmGpHmSHskmF/2lpPd1I97s9UudFLWb2qmLpHdn9ViW7bdU0qmSNulW3HXiauv3ktt/A0mLJIWkD5YZaz9xtF0PSYdIulnSGklPSLpa0mZlx9wgnnb/V8ZLukbSk9lyraQ9OxpcRHgpLMDfAS+S5jabSGqSWwO8rcn9BfwSeBh4vIp1Ae4jdQoeDBwAfJfU9/WhLus4/r0AAAWlSURBVMQ8EniINDHqAaQ57J4FvtzEvnOBe4BDSU2odwG/6OH731ZdgK8BNwDTgAmkKZCeBi6rWl0Kx5ie/V8E8MGq1YM02e5zwJey38uHga8DI6pUF2BH4CnSzCoHZcsCYCVp4GRn4uvFmzLYF9Jo/u/knm8A3A58r8n9/x5YBpwxCBJMW3UBtqlR9itgfhdiPok0bmmLXNm/AKvzZTX22zv74HpPrmyPrGz/Hr3/7dal1vs/PavLTlWqS27bkcBjwFE9TjBt/06AZ4BpvYi7w3X5FOmL54hc2cis7B86FZ+byAok7Uyat+yHfWUR8RLwI9LEmv3tvznpCrfjgBdKCrMpA6lLRNRq2vsdsEMnY6yj3qSow0iTojbab71JUUlnNP3+7krSVl0avP/Qnd9BLe3+XvqcTjqzv66E2FrRbj0+mj1+t6zA2tBuXTYG1pHOdvqsysrUqeCcYNbXN1Cz1qSZW0mqOedOzheBOyPi8o5H1rqB1qVob7ozYHW9yU0j4j7St7JaA27r7pfJT4rabe3WpZa9Sc2Ud3cmtJa1XRdJ7wCOJH3x6rV267EnqUXgKEkPSFor6beS/qq8UPvVbl0uy7Y5R9IoSaNIV/SuIH0B7QgnmPX1TS3zVKF8RWH9eiS9FTgG+GwJcbWj7boUSToSeBdp3FLZRrJ+zJDibhRzu/uVqSMxSdqOdFO+/4iIRzsUW6sGUpevA9+IiGUdj6p17dZjO9KEuKcAJwB/SzoDuFrS6E4H2aS26hIRD5H6ZA8FHsmWQ0iTAj/WqeBeE5NdShoBbN/fdhHRaLqaZvwrcHFE3D7A49TVxbrkX3Mc6QPiXyNifqeOa82R9BekZs5VwD/3OJyWSTqc9MH8t72OZYAEDAcOi4irAST9CrgX+DTpQppKkLQ96UxlMenCBUhfjq+Q9FfZWdCAvSYSDGn6mAub2E68etLM/DeDhpNmSjoQ+Gvg05K2zIo3Sau0JbAm0nQ2A1V6XV51kNSPcwWp3fxzzYc5IAOZFLVWs19/+5Wp3boA6Y+HdOXfbsBfRxqs3Cst10XSxsDZpH7JDbL/hS2y1ZtJ2jwinikj2AYG8vcVpKutAIiIlZIWA7t2MsAWtFuX40n9MB+JiLUAkq4n3drkONa/cWNbXhNNZBHxrYhQf0u2ed83/1qTZj7Z4PTxraRvN38k/WJXkE6jt8p+Pr5CdQEga5edS/qGdnhEvNiJOjShzElRu63duvQ5H5gMTO7kWWmb2qnLZsAbSE2rff8Xt2Xr5vDKhQvd1O7v5E7SF7diJ7hIfWO90G5dxgJ/6EsuABHxAvAH0s0bO+I1kWBaERF/InVkvzxppqQNsudXNdj1x6Q2zfzyXdJ15ROB/ygp5LoGUBckDQeuzJ5+MCJWlxVnDaVNitoD7dYFSSeRml4+HmkapF5rpy6rWP//4u+ydScDHysn1Iba/Z38PHuc2FeQNVmP45Wk2W3t1uVe4G1Z8ysASredfxuwvGPR9fo67sG48MrgxFNIf0wXUxicSLoEcB2wb4PjzKT342DaqgtwDeky6yOAvfJLF2IeCfwZmEeaV2466YPqy4XtlgHfLpTNBf5E6rA8mHTVT68HWrZcl+x9D+Ci4vsPbFulutQ4zhh6P9Cy3b+vy7N9/w9pcOJC0tiekVWqCykpriU1fx8EfJCUrNYC/7tj8fXiTanCQhpBvQx4HrgF2K+wfkL2TzKhwTFm0uME025dsuc1ly7FvCtplPGa7B/odGDDwjbLSRdV5Mu2zD6UnyKdPf4nNQYtdvn9b7kupC8C9X4HU6tUlxrH6GmCGeDf13Dg/wFPZPteC7y9an9fWdl+pNkinsyWhY0+z9pZPNmlmZmVwn0wZmZWCicYMzMrhROMmZmVwgnGzMxK4QRjZmalcIIxM7NSOMGYmVkpnGDMzKwU/x/m+iwhqaodRQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot estimated values for \"a\"\n",
    "plt.bar(result['values'], result['probabilities'], width=0.5/len(result['probabilities']))\n",
    "plt.xticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('\"a\" Value', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()\n",
    "\n",
    "# plot estimated values for option price (after re-scaling and reversing the c_approx-transformation)\n",
    "plt.bar(result['mapped_values'], result['probabilities'], width=1/len(result['probabilities']))\n",
    "plt.plot([exact_value, exact_value], [0,1], 'r--', linewidth=2)\n",
    "plt.xticks(size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('Estimated Option Price', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
